53988462feb01f2bc5a07fc6b153ca20fea13a5c
[u/mrichter/AliRoot.git] / HBTAN / AliHBTReaderESD.cxx
1 #include "AliHBTReaderESD.h"
2 //____________________________________________________________________
3 //////////////////////////////////////////////////////////////////////
4 //                                                                  //
5 // class AliHBTReaderESD                                            //
6 //                                                                  //
7 // reader for ALICE Event Summary Data (ESD).                       //
8 //                                                                  //
9 // Piotr.Skowronski@cern.ch                                         //
10 //                                                                  //
11 //////////////////////////////////////////////////////////////////////
12
13 #include <Riostream.h>
14
15 #include <TPDGCode.h>
16 #include <TString.h>
17 #include <TObjString.h>
18 #include <TTree.h>
19 #include <TFile.h>
20 #include <TKey.h>
21 #include <TParticle.h>
22 #include <TH1.h>
23
24 #include <AliRun.h>
25 #include <AliRunLoader.h>
26 #include <AliStack.h>
27 #include <AliESDtrack.h>
28 #include <AliKalmanTrack.h>
29 #include <AliESD.h>
30
31 #include "AliHBTRun.h"
32 #include "AliHBTEvent.h"
33 #include "AliHBTParticle.h"
34 #include "AliHBTParticleCut.h"
35 #include "AliHBTTrackPoints.h"
36 #include "AliHBTClusterMap.h"
37
38 ClassImp(AliHBTReaderESD)
39
40 AliHBTReaderESD::AliHBTReaderESD(const Char_t* esdfilename, const Char_t* galfilename):
41  fESDFileName(esdfilename),
42  fGAlFileName(galfilename),
43  fFile(0x0),
44  fRunLoader(0x0),
45  fKeyIterator(0x0),
46  fReadParticles(kFALSE),
47  fCheckParticlePID(kFALSE),
48  fNTrackPoints(0),
49  fdR(0.0),
50  fClusterMap(kFALSE),
51  fITSTrackPoints(kFALSE),
52  fMustTPC(kFALSE),
53  fNTPCClustMin(0),
54  fNTPCClustMax(150),
55  fTPCChi2PerClustMin(0.0),
56  fTPCChi2PerClustMax(10e5),
57  fChi2Min(0.0),
58  fChi2Max(10e5),
59  fC00Min(0.0),
60  fC00Max(10e5),
61  fC11Min(0.0),
62  fC11Max(10e5),
63  fC22Min(0.0),
64  fC22Max(10e5),
65  fC33Min(0.0),
66  fC33Max(10e5),
67  fC44Min(0.0),
68  fC44Max(10e5),
69  fTPCC00Min(0.0),
70  fTPCC00Max(10e5),
71  fTPCC11Min(0.0),
72  fTPCC11Max(10e5),
73  fTPCC22Min(0.0),
74  fTPCC22Max(10e5),
75  fTPCC33Min(0.0),
76  fTPCC33Max(10e5),
77  fTPCC44Min(0.0),
78  fTPCC44Max(10e5)
79
80 {
81   //cosntructor
82   if ( ((Int_t)kNSpecies) != ((Int_t)AliESDtrack::kSPECIES))
83     Fatal("AliHBTReaderESD","ESD defintions probobly changed. Ask Youra.");
84 }
85 /********************************************************************/
86   
87 AliHBTReaderESD::AliHBTReaderESD(TObjArray* dirs,const Char_t* esdfilename, const Char_t* galfilename):
88  AliHBTReader(dirs), 
89  fESDFileName(esdfilename),
90  fGAlFileName(galfilename),
91  fFile(0x0),
92  fRunLoader(0x0),
93  fKeyIterator(0x0),
94  fReadParticles(kFALSE),
95  fCheckParticlePID(kFALSE),
96  fNTrackPoints(0),
97  fdR(0.0),
98  fClusterMap(kFALSE),
99  fITSTrackPoints(kFALSE),
100  fMustTPC(kFALSE),
101  fNTPCClustMin(0),
102  fNTPCClustMax(150),
103  fTPCChi2PerClustMin(0.0),
104  fTPCChi2PerClustMax(10e5),
105  fChi2Min(0.0),
106  fChi2Max(10e5),
107  fC00Min(0.0),
108  fC00Max(10e5),
109  fC11Min(0.0),
110  fC11Max(10e5),
111  fC22Min(0.0),
112  fC22Max(10e5),
113  fC33Min(0.0),
114  fC33Max(10e5),
115  fC44Min(0.0),
116  fC44Max(10e5),
117  fTPCC00Min(0.0),
118  fTPCC00Max(10e5),
119  fTPCC11Min(0.0),
120  fTPCC11Max(10e5),
121  fTPCC22Min(0.0),
122  fTPCC22Max(10e5),
123  fTPCC33Min(0.0),
124  fTPCC33Max(10e5),
125  fTPCC44Min(0.0),
126  fTPCC44Max(10e5)
127 {
128   //cosntructor
129   if ( ((Int_t)kNSpecies) != ((Int_t)AliESDtrack::kSPECIES))
130     Fatal("AliHBTReaderESD","ESD defintions probobly changed. Ask Youra.");
131 }
132 /********************************************************************/
133
134 AliHBTReaderESD::~AliHBTReaderESD()
135 {
136  //desctructor
137   delete fRunLoader;
138   delete fKeyIterator;
139   delete fFile;
140 }
141 /**********************************************************/
142 Int_t AliHBTReaderESD::ReadNext()
143 {
144 //reads next event from fFile
145   //fRunLoader is for reading Kine
146   
147   if (AliHBTParticle::GetDebug())
148     Info("ReadNext","Entered");
149     
150   if (fParticlesEvent == 0x0)  fParticlesEvent = new AliHBTEvent();
151   if (fTracksEvent == 0x0)  fTracksEvent = new AliHBTEvent();
152   
153   fParticlesEvent->Reset();
154   fTracksEvent->Reset();
155         
156   do  //do{}while; is OK even if 0 dirs specified. In that case we try to read from "./"
157    {
158      if (fFile == 0x0)
159       {
160        fFile = OpenFile(fCurrentDir);//rl is opened here
161        if (fFile == 0x0)
162          {
163            Error("ReadNext","Cannot get fFile for dir no. %d",fCurrentDir);
164            fCurrentDir++;
165            continue;
166          }
167        fCurrentEvent = 0;
168        fKeyIterator = new TIter(fFile->GetListOfKeys());  
169 //       fFile->Dump();
170 //       fFile->GetListOfKeys()->Print();
171       } 
172      TKey* key = (TKey*)fKeyIterator->Next();
173      if (key == 0x0)
174       {
175         if (AliHBTParticle::GetDebug() > 2 )
176           {
177             Info("ReadNext","No more keys.");
178           }
179         fCurrentDir++;
180         delete fKeyIterator;
181         fKeyIterator = 0x0;
182         delete fFile;//we have to assume there is no more ESD objects in the fFile
183         fFile = 0x0;
184         delete fRunLoader;
185         fRunLoader = 0x0;
186         continue;
187       }
188      //try to read
189      
190      
191 //     TObject* esdobj = key->ReadObj();
192 //     if (esdobj == 0x0)
193 //      {
194 //        if (AliHBTParticle::GetDebug() > 2 )
195 //          {
196 //            Info("ReadNext","Key read NULL. Key Name is %s",key->GetName());
197 //            key->Dump();
198 //          }
199 //        continue;
200 //      }
201 //     esdobj->Dump();
202 //     AliESD* esd = dynamic_cast<AliESD*>(esdobj);
203      
204      TString esdname = "ESD";
205      esdname+=fCurrentEvent;
206      AliESD* esd = dynamic_cast<AliESD*>(fFile->Get(esdname));
207      if (esd == 0x0)
208       {
209 //        if (AliHBTParticle::GetDebug() > 2 )
210 //          {
211 //            Info("ReadNext","This key is not an AliESD object %s",key->GetName());
212 //          }
213         if (AliHBTParticle::GetDebug() > 2 )
214           {
215             Info("ReadNext","Can not find AliESD object named %s",esdname.Data());
216           }
217         fCurrentDir++;
218         delete fKeyIterator;
219         fKeyIterator = 0x0;
220         delete fFile;//we have to assume there is no more ESD objects in the fFile
221         fFile = 0x0;
222         delete fRunLoader;
223         fRunLoader = 0x0;
224         continue;
225       }
226      
227      ReadESD(esd);
228       
229      fCurrentEvent++;
230      fNEventsRead++;
231      delete esd;
232      return 0;//success -> read one event
233    }while(fCurrentDir < GetNumberOfDirs());//end of loop over directories specified in fDirs Obj Array  
234    
235   return 1; //no more directories to read
236 }
237 /**********************************************************/
238
239 Int_t AliHBTReaderESD::ReadESD(AliESD* esd)
240 {
241   //****** Tentative particle type "concentrations"
242   static const Double_t concentr[5]={0.05, 0., 0.85, 0.10, 0.05};
243   
244   Double_t pidtable[kNSpecies];//array used for reading pid probabilities from ESD track
245   Double_t w[kNSpecies];
246   Double_t mom[3];//momentum
247   Double_t pos[3];//position
248   Double_t vertexpos[3];//vertex position
249   //Reads one ESD
250   if (esd == 0x0)
251    {
252      Error("ReadESD","ESD is NULL");
253      return 1;
254    }
255
256   TDatabasePDG* pdgdb = TDatabasePDG::Instance();
257   if (pdgdb == 0x0)
258    {
259      Error("ReadESD","Can not get PDG Database Instance.");
260      return 1;
261    }
262    
263   Float_t mf = esd->GetMagneticField(); 
264
265   if ( (mf == 0.0) && ((fNTrackPoints > 0) || fITSTrackPoints) )
266    {
267       Error("ReadESD","Magnetic Field is 0 and Track Points Demended. Skipping to next event.");
268       return 1;
269    }
270   if (fITSTrackPoints)
271    {
272      Info("ReadESD","Magnetic Field is %f",mf/10.);
273      AliKalmanTrack::SetMagneticField(mf/10.);
274    }
275    
276   AliStack* stack = 0x0;
277   if (fReadParticles && fRunLoader)
278    {
279      fRunLoader->GetEvent(fCurrentEvent);
280      stack = fRunLoader->Stack();
281    }
282
283   const AliESDVertex* vertex = esd->GetVertex();
284   if (vertex == 0x0)
285    {
286      Info("ReadESD","ESD returned NULL pointer to vertex - assuming (0.0,0.0,0.0)");
287      vertexpos[0] = 0.0;
288      vertexpos[1] = 0.0;
289      vertexpos[2] = 0.0;
290    }
291   else
292    {
293      vertex->GetXYZ(vertexpos);
294    }
295    
296   if (AliHBTParticle::GetDebug() > 0)
297    {
298      Info("ReadESD","Primary Vertex is (%f,%f,%f)",vertexpos[0],vertexpos[1],vertexpos[2]);
299    }
300    
301   Info("ReadESD","Reading Event %d",fCurrentEvent);
302
303   Int_t ntr = esd->GetNumberOfTracks();
304   Info("ReadESD","Found %d tracks.",ntr);
305   for (Int_t i = 0;i<ntr; i++)
306    {
307      AliESDtrack *esdtrack = esd->GetTrack(i);
308      if (esdtrack == 0x0)
309       {
310         Error("Next","Can not get track %d", i);
311         continue;
312       }
313
314      //if (esdtrack->HasVertexParameters() == kFALSE) 
315      if ((esdtrack->GetStatus() & AliESDtrack::kITSrefit) == kFALSE)
316       {
317         if (AliHBTParticle::GetDebug() > 2) 
318           Info("ReadNext","Particle skipped: Data at vertex not available.");
319         continue;
320       }
321
322      if (fMustTPC)
323       {
324        if ((esdtrack->GetStatus() & AliESDtrack::kTPCin) == kFALSE)
325         {
326           if (AliHBTParticle::GetDebug() > 2) 
327             Info("ReadNext","Particle skipped: Was not reconstructed in TPC.");
328           continue;
329         }
330       }     
331      if ((esdtrack->GetStatus() & AliESDtrack::kESDpid) == kFALSE) 
332       {
333         if (AliHBTParticle::GetDebug() > 2) 
334           Info("ReadNext","Particle skipped: PID BIT is not set.");
335         continue;
336       }
337
338
339      Double_t extx;
340      Double_t extp[5];
341      esdtrack->GetConstrainedExternalParameters(extx,extp);
342      if (extp[4] == 0.0)
343       {
344         if (AliHBTParticle::GetDebug() > 2) 
345           Info("ReadNext","Track has 0 contrianed curvature -> Probobly parameters never updated. Skipping.");
346         continue;
347       } 
348      esdtrack->GetESDpid(pidtable);
349      esdtrack->GetConstrainedPxPyPz(mom);
350      esdtrack->GetConstrainedXYZ(pos);
351      pos[0] -= vertexpos[0];//we are interested only in relative position to Primary vertex at this point
352      pos[1] -= vertexpos[1];
353      pos[2] -= vertexpos[2];
354
355      Int_t charge = (extp[4] > 0)?1:-1;//if curvature=charg/Pt is positive charge is positive
356
357      //Particle from kinematics
358      AliHBTParticle* particle = 0;
359      Bool_t keeppart = kFALSE;
360      if ( fReadParticles && stack  )
361       {
362         if (esdtrack->GetLabel() < 0) continue;//this is fake -  we are not able to match any track 
363         TParticle *p = stack->Particle(esdtrack->GetLabel());
364         if (p==0x0) 
365          {
366            Error("ReadNext","Can not find track with such label.");
367            continue;
368          }
369         if (fCheckParticlePID)
370          {
371            if(Pass(p->GetPdgCode())) 
372             {
373               if ( AliHBTParticle::GetDebug() > 5 )
374                 Info("ReadNext","Simulated Particle PID (%d) did not pass the cut.",p->GetPdgCode());
375               continue; //check if we are intersted with particles of this type 
376             }
377          }
378 //           if(p->GetPdgCode()<0) charge = -1;
379         particle = new AliHBTParticle(*p,i);
380
381       }
382       
383      if(CheckTrack(esdtrack)) continue;
384       
385      //Here we apply Bayes' formula
386      Double_t rc=0.;
387      for (Int_t s=0; s<AliESDtrack::kSPECIES; s++) rc+=concentr[s]*pidtable[s];
388      if (rc==0.0) 
389       {
390         if (AliHBTParticle::GetDebug() > 2) 
391           Info("ReadNext","Particle rejected since total bayessian PID probab. is zero.");
392         continue;
393       }
394
395      for (Int_t s=0; s<AliESDtrack::kSPECIES; s++) w[s]=concentr[s]*pidtable[s]/rc;
396
397      if (AliHBTParticle::GetDebug() > 4)
398       { 
399         Info("ReadNext","###########################################################################");
400         Info("ReadNext","Momentum: %f %f %f",mom[0],mom[1],mom[2]);
401         Info("ReadNext","Position: %f %f %f",pos[0],pos[1],pos[2]);
402         TString msg("Pid list got from track:");
403         for (Int_t s = 0;s<kNSpecies;s++)
404          {
405            msg+="\n    ";
406            msg+=s;
407            msg+="(";
408            msg+=charge*GetSpeciesPdgCode((ESpecies)s);
409            msg+="): ";
410            msg+=w[s];
411            msg+=" (";
412            msg+=pidtable[s];
413            msg+=")";
414          }
415         Info("ReadNext","%s",msg.Data());
416       }//if (AliHBTParticle::GetDebug()>4)
417
418       AliHBTTrackPoints* tpts = 0x0;
419       if (fNTrackPoints > 0) 
420        {
421          tpts = new AliHBTTrackPoints(fNTrackPoints,esdtrack,mf,fdR);
422        }
423       
424       AliHBTTrackPoints* itstpts = 0x0;
425       if (fITSTrackPoints) 
426        {
427          itstpts = new AliHBTTrackPoints(AliHBTTrackPoints::kITS,esdtrack);
428        }
429       
430
431       AliHBTClusterMap* cmap = 0x0; 
432       if (  fClusterMap ) 
433        {
434          cmap = new AliHBTClusterMap(esdtrack);
435        }
436
437      for (Int_t s = 0; s<kNSpecies; s++)
438       {
439         Int_t pdgcode = charge*GetSpeciesPdgCode((ESpecies)s);
440         Float_t pp = w[s];
441         if (pp == 0.0) 
442          {
443            if ( AliHBTParticle::GetDebug() > 5 )
444              Info("ReadNext","Probability of being PID %d is zero. Continuing.",pdgcode);
445            continue;
446          }
447
448         if(Pass(pdgcode)) 
449          {
450            if ( AliHBTParticle::GetDebug() > 5 )
451              Info("ReadNext","PID (%d) did not pass the cut.",pdgcode);
452            continue; //check if we are intersted with particles of this type 
453          }
454
455         Double_t mass = pdgdb->GetParticle(pdgcode)->Mass();
456         Double_t tEtot = TMath::Sqrt( mom[0]*mom[0] + mom[1]*mom[1] + mom[2]*mom[2] + mass*mass);//total energy of the track
457
458         AliHBTParticle* track = new AliHBTParticle(pdgcode, w[s],i, 
459                                                    mom[0], mom[1], mom[2], tEtot,
460                                                    pos[0], pos[1], pos[2], 0.);
461         //copy probabilitis of other species (if not zero)
462         for (Int_t k = 0; k<kNSpecies; k++)
463          {
464            if (k == s) continue;
465            if (w[k] == 0.0) continue;
466            track->SetPIDprobability(charge*GetSpeciesPdgCode( (ESpecies)k ),w[k]);
467          }
468
469         if(Pass(track))//check if meets all criteria of any of our cuts
470                        //if it does not delete it and take next good track
471          { 
472            if ( AliHBTParticle::GetDebug() > 4 )
473              Info("ReadNext","Track did not pass the cut");
474            delete track;
475            continue;
476          }
477
478         //Single Particle cuts on cluster map and track points rather do not have sense
479         if (tpts)
480          {
481            track->SetTrackPoints(tpts); 
482          }
483
484         if (itstpts)
485          {
486            track->SetITSTrackPoints(itstpts); 
487          }
488
489         if (cmap) 
490          { 
491            track->SetClusterMap(cmap);
492          }
493
494         fTracksEvent->AddParticle(track);
495         if (particle) fParticlesEvent->AddParticle(particle);
496         keeppart = kTRUE;
497
498         if (AliHBTParticle::GetDebug() > 4 )
499          {
500            Info("ReadNext","\n\nAdding Particle with incarnation %d",pdgcode);
501            track->Print();
502            if (particle) particle->Print();
503            Info("ReadNext","\n----------------------------------------------\n");
504          }
505       }//for (Int_t s = 0; s<kNSpecies; s++)
506
507      if (keeppart == kFALSE) 
508       {
509         delete particle;//particle was not stored in event
510         delete tpts;
511         delete itstpts;
512         delete cmap;
513       }
514      else
515       {
516         if (particle->P() < 0.00001)
517          {
518            Info("ReadNext","###################################");
519            Info("ReadNext","###################################");
520            Info("ReadNext","Track Label %d",esdtrack->GetLabel());
521            TParticle *p = stack->Particle(esdtrack->GetLabel());
522            Info("ReadNext","");
523            p->Print();
524            Info("ReadNext","");
525            particle->Print();
526            TString command("touch BadInput");
527            ofstream sfile("BadEvent",ios::out);
528            sfile<<fCurrentDir<<endl;
529          }
530       } 
531
532    }//for (Int_t i = 0;i<ntr; i++)  -- loop over tracks
533
534   Info("ReadNext","Read %d tracks and %d particles from event %d (event %d in dir %d).",
535          fTracksEvent->GetNumberOfParticles(), fParticlesEvent->GetNumberOfParticles(),
536          fNEventsRead,fCurrentEvent,fCurrentDir);
537   fTrackCounter->Fill(fTracksEvent->GetNumberOfParticles());
538   return 0;
539 }
540
541 /**********************************************************/
542
543 void AliHBTReaderESD::Rewind()
544 {
545   //rewinds reading 
546   delete fKeyIterator;
547   delete fFile;
548   fFile = 0x0;
549   fKeyIterator = 0x0;
550   delete fRunLoader;
551   fRunLoader = 0x0;
552   fCurrentDir = 0;
553   fNEventsRead = 0;
554   if (fTrackCounter) fTrackCounter->Reset();
555 }
556 /**********************************************************/
557
558 TFile* AliHBTReaderESD::OpenFile(Int_t n)
559 {
560 //opens fFile with kine tree
561
562  const TString& dirname = GetDirName(n);
563  if (dirname == "")
564   {
565    Error("OpenFiles","Can not get directory name");
566    return 0x0;
567   }
568  TString filename = dirname +"/"+ fESDFileName;
569  TFile *ret = TFile::Open(filename.Data()); 
570
571  if ( ret == 0x0)
572   {
573     Error("OpenFiles","Can't open fFile %s",filename.Data());
574     return 0x0;
575   }
576  if (!ret->IsOpen())
577   {
578     Error("OpenFiles","Can't open fFile  %s",filename.Data());
579     return 0x0;
580   }
581  
582  if (fReadParticles )
583   {
584    fRunLoader = AliRunLoader::Open(dirname +"/"+ fGAlFileName);
585    if (fRunLoader == 0x0)
586     {
587       Error("OpenFiles","Can't get RunLoader for directory %s",dirname.Data());
588       delete ret;
589       return 0x0;
590     }
591     
592    fRunLoader->LoadHeader();
593    if (fRunLoader->LoadKinematics())
594     {
595       Error("Next","Error occured while loading kinematics.");
596       delete fRunLoader;
597       delete ret;
598       return 0x0;
599     }
600   }
601    
602  return ret;
603 }
604 /**********************************************************/
605
606 Int_t AliHBTReaderESD::GetSpeciesPdgCode(ESpecies spec)//skowron
607 {
608   //returns pdg code from the PID index
609   //ask jura about charge
610   switch (spec)
611    {
612      case kESDElectron:
613        return kPositron;
614        break;
615      case kESDMuon:
616        return kMuonPlus;
617        break;
618      case kESDPion:
619        return kPiPlus;
620        break;
621      case kESDKaon:
622        return kKPlus;
623        break;
624      case kESDProton:
625        return kProton;
626        break;
627      default:
628        ::Warning("GetSpeciesPdgCode","Specie with number %d is not defined.",(Int_t)spec);
629        break;
630    }
631   return 0;
632 }
633 /********************************************************************/
634 Bool_t AliHBTReaderESD::CheckTrack(AliESDtrack* t) const
635 {
636   //Performs check of the track
637   
638   if ( (t->GetConstrainedChi2() < fChi2Min) || (t->GetConstrainedChi2() > fChi2Max) ) return kTRUE;
639   
640   if ( (t->GetTPCclusters(0x0) < fNTPCClustMin) || (t->GetTPCclusters(0x0) > fNTPCClustMax) ) return kTRUE;
641
642   if (t->GetTPCclusters(0x0) > 0)
643    {
644      Float_t chisqpercl = t->GetTPCchi2()/((Double_t)t->GetTPCclusters(0x0));
645      if ( (chisqpercl < fTPCChi2PerClustMin) || (chisqpercl > fTPCChi2PerClustMax) ) return kTRUE;
646    }
647
648   Double_t cc[15];
649   t->GetConstrainedExternalCovariance(cc);
650
651   if ( (cc[0]  < fC00Min) || (cc[0]  > fC00Max) ) return kTRUE;
652   if ( (cc[2]  < fC11Min) || (cc[2]  > fC11Max) ) return kTRUE;
653   if ( (cc[5]  < fC22Min) || (cc[5]  > fC22Max) ) return kTRUE;
654   if ( (cc[9]  < fC33Min) || (cc[9]  > fC33Max) ) return kTRUE;
655   if ( (cc[14] < fC44Min) || (cc[14] > fC44Max) ) return kTRUE;
656
657
658   t->GetInnerExternalCovariance(cc);
659
660   if ( (cc[0]  < fTPCC00Min) || (cc[0]  > fTPCC00Max) ) return kTRUE;
661   if ( (cc[2]  < fTPCC11Min) || (cc[2]  > fTPCC11Max) ) return kTRUE;
662   if ( (cc[5]  < fTPCC22Min) || (cc[5]  > fTPCC22Max) ) return kTRUE;
663   if ( (cc[9]  < fTPCC33Min) || (cc[9]  > fTPCC33Max) ) return kTRUE;
664   if ( (cc[14] < fTPCC44Min) || (cc[14] > fTPCC44Max) ) return kTRUE;
665
666   return kFALSE;
667
668 }
669 /********************************************************************/
670
671 void AliHBTReaderESD::SetChi2Range(Float_t min, Float_t max)
672 {
673   //sets range of Chi2 per Cluster
674   fChi2Min = min;
675   fChi2Max = max;
676 }
677 /********************************************************************/
678
679 void AliHBTReaderESD::SetTPCNClustersRange(Int_t min,Int_t max)
680 {
681  //sets range of Number Of Clusters that tracks have to have
682  fNTPCClustMin = min;
683  fNTPCClustMax = max;
684 }
685 /********************************************************************/
686
687 void AliHBTReaderESD::SetTPCChi2PerCluserRange(Float_t min, Float_t max)
688 {
689   //sets range of Chi2 per Cluster
690   fTPCChi2PerClustMin = min;
691   fTPCChi2PerClustMax = max;
692 }
693 /********************************************************************/
694
695 void AliHBTReaderESD::SetC00Range(Float_t min, Float_t max)
696 {
697  //Sets range of C00 parameter of covariance matrix of the track
698  //it defines uncertainty of the momentum
699   fC00Min = min;
700   fC00Max = max;
701 }
702 /********************************************************************/
703
704 void AliHBTReaderESD::SetC11Range(Float_t min, Float_t max)
705 {
706  //Sets range of C11 parameter of covariance matrix of the track
707  //it defines uncertainty of the momentum
708   fC11Min = min;
709   fC11Max = max;
710 }
711 /********************************************************************/
712
713 void AliHBTReaderESD::SetC22Range(Float_t min, Float_t max)
714 {
715  //Sets range of C22 parameter of covariance matrix of the track
716  //it defines uncertainty of the momentum
717   fC22Min = min;
718   fC22Max = max;
719 }
720 /********************************************************************/
721
722 void AliHBTReaderESD::SetC33Range(Float_t min, Float_t max)
723 {
724  //Sets range of C33 parameter of covariance matrix of the track
725  //it defines uncertainty of the momentum
726   fC33Min = min;
727   fC33Max = max;
728 }
729 /********************************************************************/
730
731 void AliHBTReaderESD::SetC44Range(Float_t min, Float_t max)
732 {
733  //Sets range of C44 parameter of covariance matrix of the track
734  //it defines uncertainty of the momentum
735   fC44Min = min;
736   fC44Max = max;
737 }