Update DataLoader and ESDtrack
[u/mrichter/AliRoot.git] / STEER / AliESDtrack.cxx
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 //-----------------------------------------------------------------
17 //           Implementation of the ESD track class
18 //   ESD = Event Summary Data
19 //   This is the class to deal with during the phisical analysis of data
20 //      Origin: Iouri Belikov, CERN
21 //      e-mail: Jouri.Belikov@cern.ch
22 //-----------------------------------------------------------------
23
24 #include "TMath.h"
25
26 #include "AliESDtrack.h"
27 #include "AliKalmanTrack.h"
28
29 ClassImp(AliESDtrack)
30
31 //_______________________________________________________________________
32 AliESDtrack::AliESDtrack() : 
33 fFlags(0),
34 fLabel(0),
35 fTrackLength(0),
36 fStopVertex(0),
37 fRalpha(0),
38 fRx(0),
39 fCalpha(0),
40 fCx(0),
41 fCchi2(1e10),
42 fIalpha(0),
43 fIx(0),
44 fOalpha(0),
45 fOx(0),
46 fITSchi2(0),
47 fITSncls(0),
48 fITSsignal(0),
49 fTPCchi2(0),
50 fTPCncls(0),
51 fTPCClusterMap(159),//number of padrows
52 fTPCsignal(0),
53 fTRDchi2(0),
54 fTRDncls(0),
55 fTRDsignal(0),
56 fTOFchi2(0),
57 fTOFindex(0),
58 fTOFsignal(-1),
59 fPHOSsignal(-1),
60 fRICHsignal(-1)
61 {
62   //
63   // The default ESD constructor 
64   //
65   for (Int_t i=0; i<kSPECIES; i++) {
66     fTrackTime[i]=0.;
67     fR[i]=1.;
68     fITSr[i]=1.;
69     fTPCr[i]=1.;
70     fTRDr[i]=1.;
71     fTOFr[i]=1.;
72     fPHOSr[i]=1.;
73     fRICHr[i]=1.;
74   }
75   fPHOSr[kSPECIES]= 1.;
76   fPHOSr[kSPECIES+1]= 1.;
77   fPHOSr[kSPECIES+2]= 1.;
78   fPHOSr[kSPECIES+3]= 1.;
79
80   fPHOSpos[0]=fPHOSpos[1]=fPHOSpos[2]=0.;
81   Int_t i;
82   for (i=0; i<5; i++)  { fRp[i]=0.; fCp[i]=0.; fIp[i]=0.; fOp[i]=0.;}
83   for (i=0; i<15; i++) { fRc[i]=0.; fCc[i]=0.; fIc[i]=0.; fOc[i]=0.;   }
84   for (i=0; i<6; i++)  { fITSindex[i]=0; }
85   for (i=0; i<180; i++){ fTPCindex[i]=0; }
86   for (i=0; i<90; i++) { fTRDindex[i]=0; }
87   fTPCLabel = 0;
88   fTRDLabel = 0;
89   fITSLabel = 0;
90
91 }
92
93 //_______________________________________________________________________
94 Double_t AliESDtrack::GetMass() const {
95   // Returns the mass of the most probable particle type
96   Float_t max=0.;
97   Int_t k=-1;
98   for (Int_t i=0; i<kSPECIES; i++) {
99     if (fR[i]>max) {k=i; max=fR[i];}
100   }
101   if (k==0) return 0.00051;
102   if (k==1) return 0.10566;
103   if (k==2||k==-1) return 0.13957;
104   if (k==3) return 0.49368;
105   if (k==4) return 0.93827;
106   Warning("GetMass()","Undefined mass !");
107   return 0.13957;
108 }
109
110 //_______________________________________________________________________
111 Bool_t AliESDtrack::UpdateTrackParams(AliKalmanTrack *t, ULong_t flags) {
112   //
113   // This function updates track's running parameters 
114   //
115   SetStatus(flags);
116   fLabel=t->GetLabel();
117
118   if (t->IsStartedTimeIntegral()) {
119     SetStatus(kTIME);
120     Double_t times[10];t->GetIntegratedTimes(times); SetIntegratedTimes(times);
121     SetIntegratedLength(t->GetIntegratedLength());
122   }
123
124   fRalpha=t->GetAlpha();
125   t->GetExternalParameters(fRx,fRp);
126   t->GetExternalCovariance(fRc);
127
128   switch (flags) {
129     
130   case kITSin: case kITSout: case kITSrefit:
131     fITSncls=t->GetNumberOfClusters();
132     fITSchi2=t->GetChi2();
133     for (Int_t i=0;i<fITSncls;i++) fITSindex[i]=t->GetClusterIndex(i);
134     fITSsignal=t->GetPIDsignal();
135     fITSLabel = t->GetLabel();
136     fITSFakeRatio = t->GetFakeRatio();
137     break;
138     
139   case kTPCin: case kTPCrefit:
140     fTPCLabel = t->GetLabel();
141     fIalpha=fRalpha;
142     fIx=fRx;
143     {
144       Int_t i;
145       for (i=0; i<5; i++) fIp[i]=fRp[i];
146       for (i=0; i<15;i++) fIc[i]=fRc[i];
147     }
148   case kTPCout:
149   
150     fTPCncls=t->GetNumberOfClusters();
151     fTPCchi2=t->GetChi2();
152     
153      {//prevrow must be declared in separate namespace, otherwise compiler cries:
154       //"jump to case label crosses initialization of `Int_t prevrow'"
155        Int_t prevrow = -1;
156        //       for (Int_t i=0;i<fTPCncls;i++) 
157        for (Int_t i=0;i<160;i++) 
158         {
159           fTPCindex[i]=t->GetClusterIndex(i);
160
161           // Piotr's Cluster Map for HBT  
162           // ### please change accordingly if cluster array is changing 
163           // to "New TPC Tracking" style (with gaps in array) 
164           Int_t idx = fTPCindex[i];
165           Int_t sect = (idx&0xff000000)>>24;
166           Int_t row = (idx&0x00ff0000)>>16;
167           if (sect > 18) row +=63; //if it is outer sector, add number of inner sectors
168
169           fTPCClusterMap.SetBitNumber(row,kTRUE);
170
171           //Fill the gap between previous row and this row with 0 bits
172           //In case  ###  pleas change it as well - just set bit 0 in case there 
173           //is no associated clusters for current "i"
174           if (prevrow < 0) 
175            {
176              prevrow = row;//if previous bit was not assigned yet == this is the first one
177            }
178           else
179            { //we don't know the order (inner to outer or reverse)
180              //just to be save in case it is going to change
181              Int_t n = 0, m = 0;
182              if (prevrow < row)
183               {
184                 n = prevrow;
185                 m = row;
186               }
187              else
188               {
189                 n = row;
190                 m = prevrow;
191               }
192
193              for (Int_t j = n+1; j < m; j++)
194               {
195                 fTPCClusterMap.SetBitNumber(j,kFALSE);
196               }
197              prevrow = row; 
198            }
199           // End Of Piotr's Cluster Map for HBT
200         }
201      }
202     fTPCsignal=t->GetPIDsignal();
203     {Double_t mass=t->GetMass();    // preliminary mass setting 
204     if (mass>0.5) fR[4]=1.;         //        used by
205     else if (mass<0.4) fR[2]=1.;    // the ITS reconstruction
206     else fR[3]=1.;}
207                      //
208     break;
209
210   case kTRDout:
211     { //requested by the PHOS  ("temporary solution")
212       Double_t r=460.;
213       if (t->PropagateTo(r,30.,0.)) {  
214          fOalpha=t->GetAlpha();
215          t->GetExternalParameters(fOx,fOp);
216          t->GetExternalCovariance(fOc);
217       }
218     }
219   case kTOFin: 
220     break;
221   case kTOFout: 
222     break;
223   case kTRDin: case kTRDrefit:
224     fTRDLabel = t->GetLabel();
225
226     fTRDncls=t->GetNumberOfClusters();
227     fTRDchi2=t->GetChi2();
228     for (Int_t i=0;i<fTRDncls;i++) fTRDindex[i]=t->GetClusterIndex(i);
229     fTRDsignal=t->GetPIDsignal();
230     break;
231   case kTRDStop:
232     break;
233   default: 
234     Error("UpdateTrackParams()","Wrong flag !\n");
235     return kFALSE;
236   }
237
238   return kTRUE;
239 }
240
241 //_______________________________________________________________________
242 void 
243 AliESDtrack::SetConstrainedTrackParams(AliKalmanTrack *t, Double_t chi2) {
244   //
245   // This function sets the constrained track parameters 
246   //
247   fCalpha=t->GetAlpha();
248   t->GetExternalParameters(fCx,fCp);
249   t->GetExternalCovariance(fCc);
250   fCchi2=chi2;
251 }
252
253
254 //_______________________________________________________________________
255 void AliESDtrack::GetExternalParameters(Double_t &x, Double_t p[5]) const {
256   //---------------------------------------------------------------------
257   // This function returns external representation of the track parameters
258   //---------------------------------------------------------------------
259   x=fRx;
260   for (Int_t i=0; i<5; i++) p[i]=fRp[i];
261 }
262 //_______________________________________________________________________
263 void AliESDtrack::GetExternalCovariance(Double_t cov[15]) const {
264   //---------------------------------------------------------------------
265   // This function returns external representation of the cov. matrix
266   //---------------------------------------------------------------------
267   for (Int_t i=0; i<15; i++) cov[i]=fRc[i];
268 }
269
270
271 //_______________________________________________________________________
272 void 
273 AliESDtrack::GetConstrainedExternalParameters(Double_t &x, Double_t p[5])const{
274   //---------------------------------------------------------------------
275   // This function returns the constrained external track parameters
276   //---------------------------------------------------------------------
277   x=fCx;
278   for (Int_t i=0; i<5; i++) p[i]=fCp[i];
279 }
280 //_______________________________________________________________________
281 void 
282 AliESDtrack::GetConstrainedExternalCovariance(Double_t c[15]) const {
283   //---------------------------------------------------------------------
284   // This function returns the constrained external cov. matrix
285   //---------------------------------------------------------------------
286   for (Int_t i=0; i<15; i++) c[i]=fCc[i];
287 }
288
289
290 Double_t AliESDtrack::GetP() const {
291   //---------------------------------------------------------------------
292   // This function returns the track momentum
293   //---------------------------------------------------------------------
294   if (TMath::Abs(fRp[4])<=0) return 0;
295   Double_t pt=1./TMath::Abs(fRp[4]);
296   return pt*TMath::Sqrt(1.+ fRp[3]*fRp[3]);
297 }
298
299 void AliESDtrack::GetConstrainedPxPyPz(Double_t *p) const {
300   //---------------------------------------------------------------------
301   // This function returns the constrained global track momentum components
302   //---------------------------------------------------------------------
303   if (TMath::Abs(fCp[4])<=0) {
304     p[0]=p[1]=p[2]=0;
305     return;
306   }
307   Double_t phi=TMath::ASin(fCp[2]) + fCalpha;
308   Double_t pt=1./TMath::Abs(fCp[4]);
309   p[0]=pt*TMath::Cos(phi); p[1]=pt*TMath::Sin(phi); p[2]=pt*fCp[3]; 
310 }
311 void AliESDtrack::GetConstrainedXYZ(Double_t *xyz) const {
312   //---------------------------------------------------------------------
313   // This function returns the global track position
314   //---------------------------------------------------------------------
315   Double_t phi=TMath::ATan2(fCp[0],fCx) + fCalpha;
316   Double_t r=TMath::Sqrt(fCx*fCx + fCp[0]*fCp[0]);
317   xyz[0]=r*TMath::Cos(phi); xyz[1]=r*TMath::Sin(phi); xyz[2]=fCp[1]; 
318 }
319
320 void AliESDtrack::GetPxPyPz(Double_t *p) const {
321   //---------------------------------------------------------------------
322   // This function returns the global track momentum components
323   //---------------------------------------------------------------------
324   if (TMath::Abs(fRp[4])<=0) {
325     p[0]=p[1]=p[2]=0;
326     return;
327   }
328   Double_t phi=TMath::ASin(fRp[2]) + fRalpha;
329   Double_t pt=1./TMath::Abs(fRp[4]);
330   p[0]=pt*TMath::Cos(phi); p[1]=pt*TMath::Sin(phi); p[2]=pt*fRp[3]; 
331 }
332 void AliESDtrack::GetXYZ(Double_t *xyz) const {
333   //---------------------------------------------------------------------
334   // This function returns the global track position
335   //---------------------------------------------------------------------
336   Double_t phi=TMath::ATan2(fRp[0],fRx) + fRalpha;
337   Double_t r=TMath::Sqrt(fRx*fRx + fRp[0]*fRp[0]);
338   xyz[0]=r*TMath::Cos(phi); xyz[1]=r*TMath::Sin(phi); xyz[2]=fRp[1]; 
339 }
340
341
342 void AliESDtrack::GetInnerPxPyPz(Double_t *p) const {
343   //---------------------------------------------------------------------
344   // This function returns the global track momentum components
345   // af the entrance of the TPC
346   //---------------------------------------------------------------------
347   if (fIx==0) {p[0]=p[1]=p[2]=0.; return;}
348   Double_t phi=TMath::ASin(fIp[2]) + fIalpha;
349   Double_t pt=1./TMath::Abs(fIp[4]);
350   p[0]=pt*TMath::Cos(phi); p[1]=pt*TMath::Sin(phi); p[2]=pt*fIp[3]; 
351 }
352
353 void AliESDtrack::GetInnerXYZ(Double_t *xyz) const {
354   //---------------------------------------------------------------------
355   // This function returns the global track position
356   // af the entrance of the TPC
357   //---------------------------------------------------------------------
358   if (fIx==0) {xyz[0]=xyz[1]=xyz[2]=0.; return;}
359   Double_t phi=TMath::ATan2(fIp[0],fIx) + fIalpha;
360   Double_t r=TMath::Sqrt(fIx*fIx + fIp[0]*fIp[0]);
361   xyz[0]=r*TMath::Cos(phi); xyz[1]=r*TMath::Sin(phi); xyz[2]=fIp[1]; 
362 }
363
364 void AliESDtrack::GetInnerExternalParameters(Double_t &x, Double_t p[5]) const 
365 {
366   //skowron
367  //---------------------------------------------------------------------
368   // This function returns external representation of the track parameters at Inner Layer of TPC
369   //---------------------------------------------------------------------
370   x=fIx;
371   for (Int_t i=0; i<5; i++) p[i]=fIp[i];
372 }
373 void AliESDtrack::GetInnerExternalCovariance(Double_t cov[15]) const
374 {
375  //skowron
376  //---------------------------------------------------------------------
377  // This function returns external representation of the cov. matrix at Inner Layer of TPC
378  //---------------------------------------------------------------------
379  for (Int_t i=0; i<15; i++) cov[i]=fIc[i];
380  
381 }
382
383 void AliESDtrack::GetOuterPxPyPz(Double_t *p) const {
384   //---------------------------------------------------------------------
385   // This function returns the global track momentum components
386   // af the radius of the PHOS
387   //---------------------------------------------------------------------
388   if (fOx==0) {p[0]=p[1]=p[2]=0.; return;}
389   Double_t phi=TMath::ASin(fOp[2]) + fOalpha;
390   Double_t pt=1./TMath::Abs(fOp[4]);
391   p[0]=pt*TMath::Cos(phi); p[1]=pt*TMath::Sin(phi); p[2]=pt*fOp[3]; 
392 }
393
394 void AliESDtrack::GetOuterXYZ(Double_t *xyz) const {
395   //---------------------------------------------------------------------
396   // This function returns the global track position
397   // af the radius of the PHOS
398   //---------------------------------------------------------------------
399   if (fOx==0) {xyz[0]=xyz[1]=xyz[2]=0.; return;}
400   Double_t phi=TMath::ATan2(fOp[0],fOx) + fOalpha;
401   Double_t r=TMath::Sqrt(fOx*fOx + fOp[0]*fOp[0]);
402   xyz[0]=r*TMath::Cos(phi); xyz[1]=r*TMath::Sin(phi); xyz[2]=fOp[1]; 
403 }
404
405 //_______________________________________________________________________
406 void AliESDtrack::GetIntegratedTimes(Double_t *times) const {
407   // Returns the array with integrated times for each particle hypothesis
408   for (Int_t i=0; i<kSPECIES; i++) times[i]=fTrackTime[i];
409 }
410
411 //_______________________________________________________________________
412 void AliESDtrack::SetIntegratedTimes(const Double_t *times) {
413   // Sets the array with integrated times for each particle hypotesis
414   for (Int_t i=0; i<kSPECIES; i++) fTrackTime[i]=times[i];
415 }
416
417 //_______________________________________________________________________
418 void AliESDtrack::SetITSpid(const Double_t *p) {
419   // Sets values for the probability of each particle type (in ITS)
420   for (Int_t i=0; i<kSPECIES; i++) fITSr[i]=p[i];
421   SetStatus(AliESDtrack::kITSpid);
422 }
423
424 void AliESDtrack::SetITSChi2MIP(const Float_t *chi2mip){
425   for (Int_t i=0; i<6; i++) fITSchi2MIP[i]=chi2mip[i];
426 }
427 //_______________________________________________________________________
428 void AliESDtrack::GetITSpid(Double_t *p) const {
429   // Gets the probability of each particle type (in ITS)
430   for (Int_t i=0; i<kSPECIES; i++) p[i]=fITSr[i];
431 }
432
433 //_______________________________________________________________________
434 Int_t AliESDtrack::GetITSclusters(UInt_t *idx) const {
435   //---------------------------------------------------------------------
436   // This function returns indices of the assgined ITS clusters 
437   //---------------------------------------------------------------------
438   for (Int_t i=0; i<fITSncls; i++) idx[i]=fITSindex[i];
439   return fITSncls;
440 }
441
442 //_______________________________________________________________________
443 Int_t AliESDtrack::GetTPCclusters(Int_t *idx) const {
444   //---------------------------------------------------------------------
445   // This function returns indices of the assgined ITS clusters 
446   //---------------------------------------------------------------------
447   if (idx!=0)
448     for (Int_t i=0; i<180; i++) idx[i]=fTPCindex[i];  // MI I prefer some constant
449   return fTPCncls;
450 }
451
452 //_______________________________________________________________________
453 void AliESDtrack::SetTPCpid(const Double_t *p) {  
454   // Sets values for the probability of each particle type (in TPC)
455   for (Int_t i=0; i<kSPECIES; i++) fTPCr[i]=p[i];
456   SetStatus(AliESDtrack::kTPCpid);
457 }
458
459 //_______________________________________________________________________
460 void AliESDtrack::GetTPCpid(Double_t *p) const {
461   // Gets the probability of each particle type (in TPC)
462   for (Int_t i=0; i<kSPECIES; i++) p[i]=fTPCr[i];
463 }
464
465 //_______________________________________________________________________
466 Int_t AliESDtrack::GetTRDclusters(UInt_t *idx) const {
467   //---------------------------------------------------------------------
468   // This function returns indices of the assgined TRD clusters 
469   //---------------------------------------------------------------------
470   if (idx!=0)
471     for (Int_t i=0; i<90; i++) idx[i]=fTRDindex[i];  // MI I prefer some constant
472   return fTRDncls;
473 }
474
475 //_______________________________________________________________________
476 void AliESDtrack::SetTRDpid(const Double_t *p) {  
477   // Sets values for the probability of each particle type (in TRD)
478   for (Int_t i=0; i<kSPECIES; i++) fTRDr[i]=p[i];
479   SetStatus(AliESDtrack::kTRDpid);
480 }
481
482 //_______________________________________________________________________
483 void AliESDtrack::GetTRDpid(Double_t *p) const {
484   // Gets the probability of each particle type (in TRD)
485   for (Int_t i=0; i<kSPECIES; i++) p[i]=fTRDr[i];
486 }
487
488 //_______________________________________________________________________
489 void    AliESDtrack::SetTRDpid(Int_t iSpecies, Float_t p)
490 {
491   // Sets the probability of particle type iSpecies to p (in TRD)
492   fTRDr[iSpecies] = p;
493 }
494
495 Float_t AliESDtrack::GetTRDpid(Int_t iSpecies) const
496 {
497   // Returns the probability of particle type iSpecies (in TRD)
498   return fTRDr[iSpecies];
499 }
500
501 //_______________________________________________________________________
502 void AliESDtrack::SetTOFpid(const Double_t *p) {  
503   // Sets the probability of each particle type (in TOF)
504   for (Int_t i=0; i<kSPECIES; i++) fTOFr[i]=p[i];
505   SetStatus(AliESDtrack::kTOFpid);
506 }
507
508 //_______________________________________________________________________
509 void AliESDtrack::GetTOFpid(Double_t *p) const {
510   // Gets probabilities of each particle type (in TOF)
511   for (Int_t i=0; i<kSPECIES; i++) p[i]=fTOFr[i];
512 }
513
514
515
516 //_______________________________________________________________________
517 void AliESDtrack::SetPHOSpid(const Double_t *p) {  
518   // Sets the probability of each particle type (in PHOS)
519   for (Int_t i=0; i<kSPECIES+4; i++) fPHOSr[i]=p[i];
520   SetStatus(AliESDtrack::kPHOSpid);
521 }
522
523 //_______________________________________________________________________
524 void AliESDtrack::GetPHOSpid(Double_t *p) const {
525   // Gets probabilities of each particle type (in PHOS)
526   for (Int_t i=0; i<kSPECIES+4; i++) p[i]=fPHOSr[i];
527 }
528
529
530 //_______________________________________________________________________
531 void AliESDtrack::SetRICHpid(const Double_t *p) {  
532   // Sets the probability of each particle type (in RICH)
533   for (Int_t i=0; i<kSPECIES; i++) fRICHr[i]=p[i];
534   SetStatus(AliESDtrack::kRICHpid);
535 }
536
537 //_______________________________________________________________________
538 void AliESDtrack::GetRICHpid(Double_t *p) const {
539   // Gets probabilities of each particle type (in RICH)
540   for (Int_t i=0; i<kSPECIES; i++) p[i]=fRICHr[i];
541 }
542
543
544
545 //_______________________________________________________________________
546 void AliESDtrack::SetESDpid(const Double_t *p) {  
547   // Sets the probability of each particle type for the ESD track
548   for (Int_t i=0; i<kSPECIES; i++) fR[i]=p[i];
549   SetStatus(AliESDtrack::kESDpid);
550 }
551
552 //_______________________________________________________________________
553 void AliESDtrack::GetESDpid(Double_t *p) const {
554   // Gets probability of each particle type for the ESD track
555   for (Int_t i=0; i<kSPECIES; i++) p[i]=fR[i];
556 }
557