]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - ANALYSIS/AliReaderESD.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame_incremental
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (non-incremental) | history | HEAD
Record changes.
[u/mrichter/AliRoot.git] / ANALYSIS / AliReaderESD.cxx
... / ...
CommitLineData
1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/* $Id$ */
17
18//____________________________________________________________________
19//////////////////////////////////////////////////////////////////////
20// //
21// class AliReaderESD //
22// //
23// reader for ALICE Event Summary Data (ESD). //
24// //
25// Piotr.Skowronski@cern.ch //
26// //
27//////////////////////////////////////////////////////////////////////
28
29#include <TFile.h>
30#include <TGliteXmlEventlist.h>
31#include <TH1.h>
32#include <TPDGCode.h>
33#include <TParticle.h>
34#include <TString.h>
35
36#include "AliAOD.h"
37#include "AliAODParticle.h"
38#include "AliClusterMap.h"
39#include "AliESD.h"
40#include "AliESDtrack.h"
41#include "AliLog.h"
42#include "AliReaderESD.h"
43#include "AliRunLoader.h"
44#include "AliStack.h"
45
46ClassImp(AliReaderESD)
47
48AliReaderESD::AliReaderESD(const Char_t* esdfilename, const Char_t* galfilename):
49 fESDFileName(esdfilename),
50 fGAlFileName(galfilename),
51 fFile(0x0),
52 fRunLoader(0x0),
53 fKeyIterator(0x0),
54 fReadSim(kFALSE),
55 fCheckParticlePID(kFALSE),
56 fReadMostProbableOnly(kFALSE),
57 fNTrackPoints(0),
58 fdR(0.0),
59 fClusterMap(kFALSE),
60 fITSTrackPoints(kFALSE),
61 fITSTrackPointsType(AliTrackPoints::kITS),
62 fMustTPC(kFALSE),
63 fReadCentralBarrel(kTRUE),
64 fReadMuon(kFALSE),
65 fReadPHOS(kFALSE),
66 fNTPCClustMin(0),
67 fNTPCClustMax(1500),
68 fTPCChi2PerClustMin(0.0),
69 fTPCChi2PerClustMax(10e5),
70 fChi2Min(0.0),
71 fChi2Max(10e5),
72 fC00Min(0.0),
73 fC00Max(10e5),
74 fC11Min(0.0),
75 fC11Max(10e5),
76 fC22Min(0.0),
77 fC22Max(10e5),
78 fC33Min(0.0),
79 fC33Max(10e5),
80 fC44Min(0.0),
81 fC44Max(10e5),
82 fTPCC00Min(0.0),
83 fTPCC00Max(10e5),
84 fTPCC11Min(0.0),
85 fTPCC11Max(10e5),
86 fTPCC22Min(0.0),
87 fTPCC22Max(10e5),
88 fTPCC33Min(0.0),
89 fTPCC33Max(10e5),
90 fTPCC44Min(0.0),
91 fTPCC44Max(10e5)
92
93{
94 //cosntructor
95 if ( ((Int_t)kNSpecies) != ((Int_t)AliPID::kSPECIES))
96 Fatal("AliReaderESD","ESD defintions probobly changed. Ask Youra.");
97}
98/********************************************************************/
99
100AliReaderESD::AliReaderESD(TObjArray* dirs,const Char_t* esdfilename, const Char_t* galfilename):
101 AliReader(dirs),
102 fESDFileName(esdfilename),
103 fGAlFileName(galfilename),
104 fFile(0x0),
105 fRunLoader(0x0),
106 fKeyIterator(0x0),
107 fReadSim(kFALSE),
108 fCheckParticlePID(kFALSE),
109 fReadMostProbableOnly(kFALSE),
110 fNTrackPoints(0),
111 fdR(0.0),
112 fClusterMap(kFALSE),
113 fITSTrackPoints(kFALSE),
114 fITSTrackPointsType(AliTrackPoints::kITS),
115 fMustTPC(kFALSE),
116 fReadCentralBarrel(kTRUE),
117 fReadMuon(kFALSE),
118 fReadPHOS(kFALSE),
119 fNTPCClustMin(0),
120 fNTPCClustMax(150),
121 fTPCChi2PerClustMin(0.0),
122 fTPCChi2PerClustMax(10e5),
123 fChi2Min(0.0),
124 fChi2Max(10e5),
125 fC00Min(0.0),
126 fC00Max(10e5),
127 fC11Min(0.0),
128 fC11Max(10e5),
129 fC22Min(0.0),
130 fC22Max(10e5),
131 fC33Min(0.0),
132 fC33Max(10e5),
133 fC44Min(0.0),
134 fC44Max(10e5),
135 fTPCC00Min(0.0),
136 fTPCC00Max(10e5),
137 fTPCC11Min(0.0),
138 fTPCC11Max(10e5),
139 fTPCC22Min(0.0),
140 fTPCC22Max(10e5),
141 fTPCC33Min(0.0),
142 fTPCC33Max(10e5),
143 fTPCC44Min(0.0),
144 fTPCC44Max(10e5)
145{
146 //cosntructor
147 if ( ((Int_t)kNSpecies) != ((Int_t)AliPID::kSPECIES))
148 Fatal("AliReaderESD","ESD defintions probobly changed. Ask Youra.");
149}
150/********************************************************************/
151
152AliReaderESD::~AliReaderESD()
153{
154 //desctructor
155 delete fRunLoader;
156 delete fKeyIterator;
157 delete fFile;
158}
159
160/**********************************************************/
161Int_t AliReaderESD::ReadNext()
162{
163//reads next event from fFile
164 //fRunLoader is for reading Kine
165
166 AliDebug(1,"Entered");
167
168 if (fEventSim == 0x0) fEventSim = new AliAOD();
169 if (fEventRec == 0x0) fEventRec = new AliAOD();
170
171 fEventSim->Reset();
172 fEventRec->Reset();
173
174 do //do{}while; is OK even if 0 dirs specified. In that case we try to read from "./"
175 {
176 if (fFile == 0x0)
177 {
178 fFile = OpenFile(fCurrentDir);//rl is opened here
179 if (fFile == 0x0)
180 {
181 Error("ReadNext","Cannot get fFile for dir no. %d",fCurrentDir);
182 fCurrentDir++;
183 continue;
184 }
185 fCurrentEvent = 0;
186 }
187 TString esdname = "ESD";
188 esdname+=fCurrentEvent;
189 AliESD* esd = dynamic_cast<AliESD*>(fFile->Get(esdname));
190 if (esd == 0x0)
191 {
192 AliDebug(3,Form("Can not find AliESD object named %s",esdname.Data()));
193 fCurrentDir++;
194 delete fFile;//we have to assume there is no more ESD objects in the fFile
195 fFile = 0x0;
196 delete fRunLoader;
197 fRunLoader = 0x0;
198 continue;
199 }
200 ReadESD(esd);
201
202 fCurrentEvent++;
203 fNEventsRead++;
204 delete esd;
205 return 0;//success -> read one event
206 }while(fCurrentDir < GetNumberOfDirs());//end of loop over directories specified in fDirs Obj Array
207
208 return 1; //no more directories to read
209}
210
211/**********************************************************/
212Int_t AliReaderESD::ReadESD(AliESD* esd)
213{
214//Reads esd data
215 if (esd == 0x0)
216 {
217 Error("ReadESD","ESD is NULL");
218 return 1;
219 }
220
221 // seperate each method
222 if (fReadCentralBarrel) ReadESDCentral(esd);
223
224 if (fReadMuon) ReadESDMuon(esd);
225
226 if (fReadPHOS) ReadESDPHOS(esd);
227
228 return 1;
229}
230
231/**********************************************************/
232Int_t AliReaderESD::ReadESDCentral(AliESD* esd)
233{
234 //****** Tentative particle type "concentrations"
235 static const Double_t kConcentr[5]={0.05, 0., 0.85, 0.10, 0.05};
236
237 Double_t pidtable[kNSpecies];//array used for reading pid probabilities from ESD track
238 Double_t w[kNSpecies];
239 Double_t mom[3];//momentum
240 Double_t pos[3];//position
241 Double_t vertexpos[3];//vertex position
242 //Reads one ESD
243
244 TDatabasePDG* pdgdb = TDatabasePDG::Instance();
245 if (pdgdb == 0x0)
246 {
247 Error("ReadESD","Can not get PDG Database Instance.");
248 return 1;
249 }
250
251 Float_t mf = esd->GetMagneticField()/10.; //AliESD::GetMagnField returns mf in kG
252
253 if ( (mf == 0.0) && ((fNTrackPoints > 0) || fITSTrackPoints) )
254 {
255 Error("ReadESD","Magnetic Field is 0 and Track Points Demended. Skipping to next event.");
256
257 }
258
259 if (fITSTrackPoints)
260 {
261 Info("ReadESD","Magnetic Field is %f",mf);
262 //AliKalmanTrack::SetMagneticField(mf);
263 }
264
265 AliStack* stack = 0x0;
266 if (fReadSim && fRunLoader)
267 {
268 fRunLoader->GetEvent(fCurrentEvent);
269 stack = fRunLoader->Stack();
270 }
271
272 const AliESDVertex* vertex = esd->GetVertex();
273 if (vertex == 0x0)
274 {
275 Info("ReadESD","ESD returned NULL pointer to vertex - assuming (0.0,0.0,0.0)");
276 vertexpos[0] = 0.0;
277 vertexpos[1] = 0.0;
278 vertexpos[2] = 0.0;
279 }
280 else
281 {
282 vertex->GetXYZ(vertexpos);
283 }
284
285 AliDebug(1,Form("Primary Vertex is (%f,%f,%f)",vertexpos[0],vertexpos[1],vertexpos[2]));
286
287 Info("ReadESD","Reading Event %d",fCurrentEvent);
288
289 Int_t ntr = esd->GetNumberOfTracks();
290 Info("ReadESD","Found %d tracks.",ntr);
291 for (Int_t i = 0;i<ntr; i++)
292 {
293 AliESDtrack *esdtrack = esd->GetTrack(i);
294 if (esdtrack == 0x0)
295 {
296 Error("Next","Can not get track %d", i);
297 continue;
298 }
299
300 //if (esdtrack->HasVertexParameters() == kFALSE)
301 if ((esdtrack->GetStatus() & AliESDtrack::kITSrefit) == kFALSE)
302 {
303 AliDebug(3,Form("Particle skipped: Data at vertex not available."));
304 continue;
305 }
306
307 if (fMustTPC)
308 {
309 if ((esdtrack->GetStatus() & AliESDtrack::kTPCin) == kFALSE)
310 {
311 AliDebug(3,"Particle skipped: Was not reconstructed in TPC.");
312 continue;
313 }
314 }
315
316 if ((esdtrack->GetStatus() & AliESDtrack::kESDpid) == kFALSE)
317 {
318 AliDebug(3,"Particle skipped: PID BIT is not set.");
319 continue;
320 }
321
322
323 Double_t alpha,extx;
324 Double_t extp[5];
325 esdtrack->GetConstrainedExternalParameters(alpha,extx,extp);
326 if (extp[4] == 0.0)
327 {
328 AliDebug(3,"Track has 0 contrianed curvature -> Probobly parameters never updated. Skipping.");
329 continue;
330 }
331 esdtrack->GetESDpid(pidtable);
332 esdtrack->GetConstrainedPxPyPz(mom);
333 esdtrack->GetConstrainedXYZ(pos);
334 pos[0] -= vertexpos[0];//we are interested only in relative position to Primary vertex at this point
335 pos[1] -= vertexpos[1];
336 pos[2] -= vertexpos[2];
337
338 Int_t charge = (extp[4] > 0)?1:-1;//if curvature=charg/Pt is positive charge is positive
339
340 //Particle from kinematics
341 AliAODParticle* particle = 0;
342 Bool_t keeppart = kFALSE;
343 if ( fReadSim && stack )
344 {
345 if (esdtrack->GetLabel() < 0) continue;//this is fake - we are not able to match any track
346 TParticle *p = stack->Particle(esdtrack->GetLabel());
347 if (p==0x0)
348 {
349 Error("ReadNext","Can not find track with such label.");
350 continue;
351 }
352
353 if (fCheckParticlePID)
354 {
355 if(Rejected(p->GetPdgCode()))
356 {
357 AliDebug(6,Form("Simulated Particle PID (%d) did not pass the cut.",p->GetPdgCode()));
358 continue; //check if we are intersted with particles of this type
359 }
360 }
361// if(p->GetPdgCode()<0) charge = -1;
362 particle = new AliAODParticle(*p,i);
363
364 }
365
366 if(CheckTrack(esdtrack)) continue;
367
368 //Here we apply Bayes' formula
369 Double_t rc=0.;
370 for (Int_t s=0; s<AliPID::kSPECIES; s++) rc+=kConcentr[s]*pidtable[s];
371 if (rc==0.0)
372 {
373 AliDebug(3,"Particle rejected since total bayessian PID probab. is zero.");
374 continue;
375 }
376
377 for (Int_t s=0; s<AliPID::kSPECIES; s++) w[s]=kConcentr[s]*pidtable[s]/rc;
378
379 if (AliDebugLevel() > 4)
380 {
381 Info("ReadNext","###########################################################################");
382 Info("ReadNext","Momentum: %f %f %f",mom[0],mom[1],mom[2]);
383 Info("ReadNext","Position: %f %f %f",pos[0],pos[1],pos[2]);
384 TString msg("Pid list got from track:");
385 for (Int_t s = 0;s<kNSpecies;s++)
386 {
387 msg+="\n ";
388 msg+=s;
389 msg+="(";
390 msg+=charge*GetSpeciesPdgCode((ESpecies)s);
391 msg+="): ";
392 msg+=w[s];
393 msg+=" (";
394 msg+=pidtable[s];
395 msg+=")";
396 }
397 Info("ReadNext","%s",msg.Data());
398 }//if (AliDebugLevel()>4)
399
400 AliTrackPoints* tpts = 0x0;
401 if (fNTrackPoints > 0)
402 {
403 tpts = new AliTrackPoints(fNTrackPoints,esdtrack,mf*10.0,fdR);
404// tpts->Move(-vertexpos[0],-vertexpos[1],-vertexpos[2]);
405 }
406
407 AliTrackPoints* itstpts = 0x0;
408 if (fITSTrackPoints)
409 {
410 itstpts = new AliTrackPoints(fITSTrackPointsType,esdtrack,mf*10.0);
411// itstpts->Move(-vertexpos[0],-vertexpos[1],-vertexpos[2]);
412 }
413
414 AliClusterMap* cmap = 0x0;
415 if ( fClusterMap )
416 {
417 cmap = new AliClusterMap(esdtrack);
418 }
419
420 //If this flag fReadMostProbableOnly is false the
421 //loop over species (see "LOOPPIDS") is over all possible PIDs
422 //in other case the most probablle PID is searched
423 //and the loop is limited to that PID only
424
425 Int_t firstspecie = 0;
426 Int_t lastspecie = kNSpecies;
427
428 if (fReadMostProbableOnly)
429 {
430 //find the most probable PID
431 Int_t spec = 0;
432 Float_t maxprob = w[0];
433 for (Int_t s=1; s<AliPID::kSPECIES; s++)
434 {
435 if (w[s]>maxprob)
436 {
437 maxprob = w[s];
438 spec = s;
439 }
440 }
441 firstspecie = spec;
442 lastspecie = spec + 1;
443 }
444
445 for (Int_t s = firstspecie; s<lastspecie; s++)//LOOPPIDS
446 {
447 Int_t pdgcode = charge*GetSpeciesPdgCode((ESpecies)s);
448 Float_t pp = w[s];
449 if (pp == 0.0)
450 {
451 AliDebug(6,Form("Probability of being PID %d is zero. Continuing.",pdgcode));
452 continue;
453 }
454
455 if(Rejected(pdgcode))
456 {
457 AliDebug(6,Form("PID (%d) did not pass the cut.",pdgcode));
458 continue; //check if we are intersted with particles of this type
459 }
460
461 Double_t mass = pdgdb->GetParticle(pdgcode)->Mass();
462 Double_t tEtot = TMath::Sqrt( mom[0]*mom[0] + mom[1]*mom[1] + mom[2]*mom[2] + mass*mass);//total energy of the track
463
464 AliAODParticle* track = new AliAODParticle(pdgcode, w[s],i,
465 mom[0], mom[1], mom[2], tEtot,
466 pos[0], pos[1], pos[2], 0.);
467 //copy probabilitis of other species (if not zero)
468 for (Int_t k = 0; k<kNSpecies; k++)
469 {
470 if (k == s) continue;
471 if (w[k] == 0.0) continue;
472 track->SetPIDprobability(charge*GetSpeciesPdgCode( (ESpecies)k ),w[k]);
473 }
474
475 if(Rejected(track))//check if meets all criteria of any of our cuts
476 //if it does not delete it and take next good track
477 {
478 AliDebug(5,"Track did not pass the cut");
479 delete track;
480 continue;
481 }
482
483 //Single Particle cuts on cluster map and track points rather do not have sense
484 if (tpts)
485 {
486 track->SetTPCTrackPoints(tpts);
487 }
488
489 if (itstpts)
490 {
491 track->SetITSTrackPoints(itstpts);
492 }
493
494 if (cmap)
495 {
496 track->SetClusterMap(cmap);
497 }
498
499 fEventRec->AddParticle(track);
500 if (particle) fEventSim->AddParticle(particle);
501 keeppart = kTRUE;
502
503 if (AliDebugLevel() > 4 )
504 {
505 Info("ReadNext","\n\nAdding Particle with incarnation %d",pdgcode);
506 track->Print();
507 if (particle) particle->Print();
508 Info("ReadNext","\n----------------------------------------------\n");
509 }
510 }//for (Int_t s = 0; s<kNSpecies; s++)
511
512 if (keeppart == kFALSE)
513 {
514 delete particle;//particle was not stored in event
515 delete tpts;
516 delete itstpts;
517 delete cmap;
518 }
519 else
520 {
521 if ( fReadSim && stack )
522 {
523 if (particle->P() < 0.00001)
524 {
525 Info("ReadNext","###################################");
526 Info("ReadNext","###################################");
527 Info("ReadNext","Track Label %d",esdtrack->GetLabel());
528 TParticle *p = stack->Particle(esdtrack->GetLabel());
529 Info("ReadNext","");
530 p->Print();
531 Info("ReadNext","");
532 particle->Print();
533 }
534 }
535 }
536
537 }//for (Int_t i = 0;i<ntr; i++) -- loop over tracks
538
539 Info("ReadNext","Read %d tracks and %d particles from event %d (event %d in dir %d).",
540 fEventRec->GetNumberOfParticles(), fEventSim->GetNumberOfParticles(),
541 fNEventsRead,fCurrentEvent,fCurrentDir);
542 fTrackCounter->Fill(fEventRec->GetNumberOfParticles());
543
544 /******************************************************/
545 /****** Setting glevet properties *************/
546 /******************************************************/
547
548 if (fEventRec->GetNumberOfParticles() > 0)
549 {
550 fEventRec->SetPrimaryVertex(vertexpos[0],vertexpos[1],vertexpos[2]);
551 }
552 return 0;
553}
554
555/**********************************************************/
556Int_t AliReaderESD::ReadESDMuon(AliESD* esd)
557{
558 // Reads the muon tracks from the ESD
559 Double_t vertexpos[3];//vertex position, assuming no secondary decay
560
561 const AliESDVertex* vertex = esd->GetVertex();
562
563 if (vertex == 0x0) {
564 Info("ReadESD","ESD returned NULL pointer to vertex - assuming (0.0,0.0,0.0)");
565 vertexpos[0] = 0.0;
566 vertexpos[1] = 0.0;
567 vertexpos[2] = 0.0;
568 } else {
569 vertex->GetXYZ(vertexpos);
570 }
571
572 Int_t nTracks = (Int_t)esd->GetNumberOfMuonTracks() ;
573
574 AliDebug(1,Form("Reading Event %d \nFound %d tracks.",fCurrentEvent,nTracks));
575
576 // settings
577 Float_t chi2Cut = 100.;
578 Float_t ptCutMin = 1.;
579 Float_t ptCutMax = 10000.;
580 Float_t muonMass = 0.105658389;
581 Int_t pdgcode = -13;
582 Double_t thetaX, thetaY, pYZ;
583 Double_t pxRec1, pyRec1, pzRec1, e1;
584 Int_t charge;
585
586 Int_t ntrackhits;
587 Double_t fitfmin;
588
589 TLorentzVector fV1;
590 fEventRec->Reset();
591 for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
592
593 AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
594
595 thetaX = muonTrack->GetThetaX();
596 thetaY = muonTrack->GetThetaY();
597
598 pYZ = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
599 pzRec1 = - pYZ / TMath::Sqrt(1.0 + TMath::Tan(thetaY)*TMath::Tan(thetaX));
600 pxRec1 = pzRec1 * TMath::Tan(thetaX);
601 pyRec1 = pzRec1 * TMath::Tan(thetaY);
602 charge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
603 e1 = TMath::Sqrt(muonMass * muonMass + pxRec1 * pxRec1 + pyRec1 * pyRec1 + pzRec1 * pzRec1);
604 fV1.SetPxPyPzE(pxRec1, pyRec1, pzRec1, e1);
605
606 ntrackhits = muonTrack->GetNHit();
607 fitfmin = muonTrack->GetChi2();
608
609 // transverse momentum
610 Float_t pt1 = fV1.Pt();
611
612 // chi2 per d.o.f.
613 Float_t ch1 = fitfmin / (2.0 * ntrackhits - 5);
614
615 if ((ch1 < chi2Cut) && (pt1 > ptCutMin) && (pt1 < ptCutMax)) {
616 AliAODParticle* track = new AliAODParticle(pdgcode*charge,1,iTrack,
617 pxRec1, pyRec1,pzRec1, e1,
618 vertexpos[0], vertexpos[1], vertexpos[2], 0.);
619 fEventRec->AddParticle(track);
620 }
621
622 }
623 fTrackCounter->Fill(fEventRec->GetNumberOfParticles());
624 return 0;
625}
626
627/**********************************************************/
628
629void AliReaderESD::Rewind()
630{
631 //rewinds reading
632 // delete fKeyIterator;
633 delete fFile;
634 fFile = 0x0;
635 // fKeyIterator = 0x0;
636 delete fRunLoader;
637 fRunLoader = 0x0;
638 fCurrentDir = 0;
639 fNEventsRead = 0;
640 if (fEventList) fEventList->Reset();
641 if (fTrackCounter) fTrackCounter->Reset();
642}
643/**********************************************************/
644
645TFile* AliReaderESD::OpenFile(Int_t n)
646{
647//opens fFile with tree
648 if (fEventList)
649 {
650 if (fCurrentDir > n)
651 {
652 fEventList->Reset();
653 fCurrentDir = 0;
654 }
655
656 while (fCurrentDir < n)
657 {
658 fEventList->Next();
659 fCurrentDir++;
660 }
661 fEventList->Next();
662 }
663
664 const TString& dirname = GetDirName(n);
665 if (dirname == "")
666 {
667 Error("OpenFiles","Can not get directory name");
668 return 0x0;
669 }
670 TString filename;
671
672 if (fEventList)
673 {
674 filename = fEventList->GetURL(fESDFileName);
675 }
676 else
677 {
678 filename = dirname +"/"+ fESDFileName;
679 }
680 Info("OpenFile","%s ==> %s",fESDFileName.Data(),filename.Data());
681 TFile *ret = TFile::Open(filename.Data());
682
683 if ( ret == 0x0)
684 {
685 Error("OpenFiles","Can't open fFile %s",filename.Data());
686 return 0x0;
687 }
688 if (!ret->IsOpen())
689 {
690 Error("OpenFiles","Can't open fFile %s",filename.Data());
691 return 0x0;
692 }
693
694 if (fReadSim )
695 {
696 TString gafilename;
697 if (fEventList)
698 {
699 gafilename = fEventList->GetURL(fGAlFileName);
700 }
701 else
702 {
703 gafilename = dirname +"/"+ fGAlFileName;
704 }
705 Info("OpenFile","%s ==> %s",fGAlFileName.Data(),gafilename.Data());
706
707 fRunLoader = AliRunLoader::Open(gafilename);
708
709 if (fRunLoader == 0x0)
710 {
711 Error("OpenFiles","Can't get RunLoader for directory %s",dirname.Data());
712 delete ret;
713 return 0x0;
714 }
715
716 fRunLoader->LoadHeader();
717
718 if (fEventList)
719 {
720 TString kinefilename = fEventList->GetURL("Kinematics.root");
721 fRunLoader->SetKineFileName(kinefilename);
722 Info("OpenFile","%s ==> %s","Kinematics.root",kinefilename.Data());
723 }
724
725 if (fRunLoader->LoadKinematics())
726 {
727 Error("Next","Error occured while loading kinematics.");
728 delete fRunLoader;
729 delete ret;
730 return 0x0;
731 }
732 }
733
734 return ret;
735}
736/**********************************************************/
737
738Int_t AliReaderESD::GetSpeciesPdgCode(ESpecies spec)//skowron
739{
740 //returns pdg code from the PID index
741 //ask jura about charge
742 switch (spec)
743 {
744 case kESDElectron:
745 return kPositron;
746 break;
747 case kESDMuon:
748 return kMuonPlus;
749 break;
750 case kESDPion:
751 return kPiPlus;
752 break;
753 case kESDKaon:
754 return kKPlus;
755 break;
756 case kESDProton:
757 return kProton;
758 break;
759 default:
760 ::Warning("GetSpeciesPdgCode","Specie with number %d is not defined.",(Int_t)spec);
761 break;
762 }
763 return 0;
764}
765/********************************************************************/
766Bool_t AliReaderESD::CheckTrack(AliESDtrack* t) const
767{
768 //Performs check of the track
769
770 if ( (t->GetConstrainedChi2() < fChi2Min) || (t->GetConstrainedChi2() > fChi2Max) ) return kTRUE;
771
772 if ( (t->GetTPCclusters(0x0) < fNTPCClustMin) || (t->GetTPCclusters(0x0) > fNTPCClustMax) ) return kTRUE;
773
774 if (t->GetTPCclusters(0x0) > 0)
775 {
776 Float_t chisqpercl = t->GetTPCchi2()/((Double_t)t->GetTPCclusters(0x0));
777 if ( (chisqpercl < fTPCChi2PerClustMin) || (chisqpercl > fTPCChi2PerClustMax) ) return kTRUE;
778 }
779
780 Double_t cc[15];
781 t->GetConstrainedExternalCovariance(cc);
782
783 if ( (cc[0] < fC00Min) || (cc[0] > fC00Max) ) return kTRUE;
784 if ( (cc[2] < fC11Min) || (cc[2] > fC11Max) ) return kTRUE;
785 if ( (cc[5] < fC22Min) || (cc[5] > fC22Max) ) return kTRUE;
786 if ( (cc[9] < fC33Min) || (cc[9] > fC33Max) ) return kTRUE;
787 if ( (cc[14] < fC44Min) || (cc[14] > fC44Max) ) return kTRUE;
788
789
790 t->GetInnerExternalCovariance(cc);
791
792 if ( (cc[0] < fTPCC00Min) || (cc[0] > fTPCC00Max) ) return kTRUE;
793 if ( (cc[2] < fTPCC11Min) || (cc[2] > fTPCC11Max) ) return kTRUE;
794 if ( (cc[5] < fTPCC22Min) || (cc[5] > fTPCC22Max) ) return kTRUE;
795 if ( (cc[9] < fTPCC33Min) || (cc[9] > fTPCC33Max) ) return kTRUE;
796 if ( (cc[14] < fTPCC44Min) || (cc[14] > fTPCC44Max) ) return kTRUE;
797
798 return kFALSE;
799
800}
801/********************************************************************/
802
803void AliReaderESD::SetChi2Range(Float_t min, Float_t max)
804{
805 //sets range of Chi2 per Cluster
806 fChi2Min = min;
807 fChi2Max = max;
808}
809/********************************************************************/
810
811void AliReaderESD::SetTPCNClustersRange(Int_t min,Int_t max)
812{
813 //sets range of Number Of Clusters that tracks have to have
814 fNTPCClustMin = min;
815 fNTPCClustMax = max;
816}
817/********************************************************************/
818
819void AliReaderESD::SetTPCChi2PerCluserRange(Float_t min, Float_t max)
820{
821 //sets range of Chi2 per Cluster
822 fTPCChi2PerClustMin = min;
823 fTPCChi2PerClustMax = max;
824}
825/********************************************************************/
826
827void AliReaderESD::SetC00Range(Float_t min, Float_t max)
828{
829 //Sets range of C00 parameter of covariance matrix of the track
830 //it defines uncertainty of the momentum
831 fC00Min = min;
832 fC00Max = max;
833}
834/********************************************************************/
835
836void AliReaderESD::SetC11Range(Float_t min, Float_t max)
837{
838 //Sets range of C11 parameter of covariance matrix of the track
839 //it defines uncertainty of the momentum
840 fC11Min = min;
841 fC11Max = max;
842}
843/********************************************************************/
844
845void AliReaderESD::SetC22Range(Float_t min, Float_t max)
846{
847 //Sets range of C22 parameter of covariance matrix of the track
848 //it defines uncertainty of the momentum
849 fC22Min = min;
850 fC22Max = max;
851}
852/********************************************************************/
853
854void AliReaderESD::SetC33Range(Float_t min, Float_t max)
855{
856 //Sets range of C33 parameter of covariance matrix of the track
857 //it defines uncertainty of the momentum
858 fC33Min = min;
859 fC33Max = max;
860}
861/********************************************************************/
862
863void AliReaderESD::SetC44Range(Float_t min, Float_t max)
864{
865 //Sets range of C44 parameter of covariance matrix of the track
866 //it defines uncertainty of the momentum
867 fC44Min = min;
868 fC44Max = max;
869}