TOFtrackerV2 cpu usage optimezed
[u/mrichter/AliRoot.git] / TOF / AliTOFtrack.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 /* $Id$ */
17
18 /////////////////////////////////////////////////////////////////////////////
19 //                                                                         //
20 // AliTOFtrack class                                                       //
21 //                                                                         //
22 // Authors: Bologna-CERN-ITEP-Salerno Group                                //
23 //                                                                         //
24 // Description: class for handling ESD extracted tracks for TOF matching.  //
25 //                                                                         //
26 /////////////////////////////////////////////////////////////////////////////
27
28 #include "AliESDtrack.h" 
29 #include "AliTracker.h" 
30
31 #include "AliTOFGeometry.h"
32 #include "AliTOFtrack.h" 
33
34 ClassImp(AliTOFtrack)
35
36 //_____________________________________________________________________________
37 AliTOFtrack::AliTOFtrack() : 
38   AliKalmanTrack(),
39   fSeedInd(-1),
40   fSeedLab(-1)
41 {
42   //
43   // Default constructor.
44   //
45 }                                
46
47 //_____________________________________________________________________________
48 AliTOFtrack::AliTOFtrack(const AliTOFtrack& t) : 
49   AliKalmanTrack(t),
50   fSeedInd(t.fSeedInd),
51   fSeedLab(t.fSeedLab) 
52 {
53   //
54   // Copy constructor.
55   //
56 }                                
57
58 //____________________________________________________________________________
59 AliTOFtrack& AliTOFtrack::operator=(const AliESDtrack& t)
60 {
61   // ass. op.
62   SetLabel(t.GetLabel());
63   SetChi2(0.);
64   SetMass(t.GetMassForTracking());
65
66   Set(t.GetX(),t.GetAlpha(),t.GetParameter(),t.GetCovariance());
67
68   if ((t.GetStatus()&AliESDtrack::kTIME) == 0) return *this;
69   StartTimeIntegral();
70   Double_t times[10]; 
71   for(Int_t isp=0;isp<AliPID::kSPECIESC;isp++){
72     times[isp] = t.GetIntegratedTimesOld(isp); // in ps
73   }
74
75   SetIntegratedTimes(times);
76   SetIntegratedLength(t.GetIntegratedLengthOld());
77
78   return *this;
79
80 }
81
82 //_____________________________________________________________________________
83 AliTOFtrack::AliTOFtrack(const AliESDtrack& t) :
84   AliKalmanTrack(), 
85   fSeedInd(-1),
86   fSeedLab(-1) 
87 {
88   //
89   // Constructor from AliESDtrack
90   //
91   SetLabel(t.GetLabel());
92   SetChi2(0.);
93   SetMass(t.GetMassForTracking());
94
95   Set(t.GetX(),t.GetAlpha(),t.GetParameter(),t.GetCovariance());
96
97   if ((t.GetStatus()&AliESDtrack::kTIME) == 0) return;
98   StartTimeIntegral();
99   Double_t times[10]; 
100   for(Int_t isp=0;isp<AliPID::kSPECIESC;isp++){
101     times[isp] = t.GetIntegratedTimesOld(isp); // in ps
102   }
103
104   SetIntegratedTimes(times);
105   SetIntegratedLength(t.GetIntegratedLengthOld());
106
107 }              
108
109 //____________________________________________________________________________
110 AliTOFtrack& AliTOFtrack::operator=(const AliTOFtrack &/*source*/)
111 {
112   // ass. op.
113
114   return *this;
115
116 }
117
118 //_____________________________________________________________________________
119 Bool_t AliTOFtrack::PropagateTo(Double_t xk,Double_t /*x0*/,Double_t /*rho*/)
120  {
121   //
122   // Propagates a track of particle with mass=pm to a reference plane 
123   // defined by x=xk through media of density=rho and radiationLength=x0
124   //
125
126   if (xk == GetX()) return kTRUE;
127   
128   Double_t oldX=GetX();//, oldY=GetY(), oldZ=GetZ();
129   Double_t start[3], end[3], mparam[7];
130
131   /* get start position */
132   GetXYZ(start);
133
134   /* propagate the track */
135   Double_t b[3];GetBxByBz(b);
136   if (!AliExternalTrackParam::PropagateToBxByBz(xk,b)) return kFALSE;
137   // OLD used code
138   //Double_t bz=GetBz();
139   //if (!AliExternalTrackParam::PropagateTo(xk,bz)) return kFALSE;
140
141   /* get end position */
142   GetXYZ(end);
143
144   /* add time step to integral */
145 #if 0 /*** OLD ***/
146   Double_t d = TMath::Sqrt((GetX()-oldX)*(GetX()-oldX) + 
147                            (GetY()-oldY)*(GetY()-oldY) + 
148                            (GetZ()-oldZ)*(GetZ()-oldZ));
149   if (IsStartedTimeIntegral() && GetX()>oldX) AddTimeStep(d);
150 #endif
151   Double_t d = TMath::Sqrt((end[0]-start[0])*(end[0]-start[0]) + 
152                            (end[1]-start[1])*(end[1]-start[1]) + 
153                            (end[2]-start[2])*(end[2]-start[2]));
154   if (IsStartedTimeIntegral() && GetX()>oldX) AddTimeStep(d);
155
156   /* get material budget from tracker */
157   AliTracker::MeanMaterialBudget(start, end, mparam);
158   Double_t xTimesRho = mparam[4]*mparam[0];
159   if (oldX < xk) { // CZ
160     xTimesRho = -xTimesRho; // it should be negative in case of outward
161                             // propagation (--> energy decreases)
162   } // CZ
163   Double_t xOverX0   = mparam[1];
164
165   /* correct for mean material */
166   if (!AliExternalTrackParam::CorrectForMeanMaterial(xOverX0,xTimesRho,GetMass())) return kFALSE;
167
168
169 #if 0 /*** OLD ***/
170   if (!AliExternalTrackParam::CorrectForMaterial(d*rho/x0,x0,GetMass())) 
171      return kFALSE;
172 #endif
173
174   /*
175   //Energy losses************************
176   if((5940*beta2/(1-beta2+1e-10) - beta2) < 0){return 0;}
177
178   Double_t dE=0.153e-3/beta2*(log(5940*beta2/(1-beta2+1e-10)) - beta2)*d*rho;
179   //
180   // suspicious part - think about it ?
181   Double_t kinE =  TMath::Sqrt(p2);
182   if (dE>0.8*kinE) dE = 0.8*kinE;  //      
183   if (dE<0)        dE = 0.0;       // not valid region for Bethe bloch 
184   */
185
186   return kTRUE;            
187 }     
188
189 //_____________________________________________________________________________
190 Bool_t AliTOFtrack::PropagateToInnerTOF()
191 {
192   // Propagates a track of particle with mass=pm to a reference plane 
193   // defined by x=xk through media of density=rho and radiationLength=x0
194
195   //const Double_t kAlphac  = TMath::Pi()/9.0; // 20 degree
196   const Double_t kAlphac  = AliTOFGeometry::GetAlpha(); // 20 degree
197   const Double_t kTalphac = TMath::Tan(kAlphac*0.5);
198
199   //const Double_t kStepSize = 0.1; // [cm] Step size
200   const Double_t kStepSize = 0.5; // [cm] Step size
201
202   Double_t x = GetX();
203   //Double_t bz = GetBz();
204
205   //Double_t xyz0[3];
206   //Double_t xyz1[3];
207   //Double_t y;
208   //Double_t z;
209
210   Int_t nsteps = (Int_t)((AliTOFGeometry::Rmin()-x)/kStepSize);
211   for (Int_t istep=0; istep<nsteps; istep++){
212
213     // Critical alpha  - cross sector indication
214
215     Double_t dir = (GetX() > AliTOFGeometry::Rmin()) ? -1.0 : 1.0;
216
217     x = GetX()+dir*kStepSize;
218     if ( x<GetX() && GetX()<AliTOFGeometry::Rmin()) {
219       AliDebug(1,Form("Track doesn't reach rho=%f",AliTOFGeometry::Rmin()));
220       return kFALSE;
221     }
222
223     //GetXYZ(xyz0);
224     //bz = GetBz();
225     //GetXYZAt(x,bz,xyz1);
226     //AliExternalTrackParam::GetYAt(x,bz,y);
227     //AliExternalTrackParam::GetZAt(x,bz,z);
228  
229     if (!(PropagateTo(x,0.,0.))) { /* passing 0.,0. as arguments since now
230                                       this method queries TGeo for material budget
231                                    */
232       return kFALSE;
233     }
234     
235     if (GetY() >  GetX()*kTalphac) {
236       if (!(Rotate(kAlphac))) return kFALSE;
237     } else if (GetY() < -GetX()*kTalphac) {
238       if (!(Rotate(-kAlphac))) return kFALSE;
239     }
240
241   }
242
243   //Bool_t check = PropagateTo(AliTOFGeometry::RinTOF());
244   Bool_t check = PropagateTo(AliTOFGeometry::RinTOF(),0.,0.); /* passing 0.,0. as arguments since now
245                                                                  this method queries TGeo for material budget
246                                                               */
247
248   if (!check) return kFALSE;
249
250   if (GetY() >  GetX()*kTalphac) {
251     if (!(Rotate(kAlphac))) return kFALSE;
252   } else if (GetY() < -GetX()*kTalphac) {
253     if (!(Rotate(-kAlphac))) return kFALSE;
254   }
255
256   return kTRUE;
257   
258 }     
259
260 //_____________________________________________________________________________
261 Bool_t AliTOFtrack::PropagateToInnerTOFold()
262 {
263   // Propagates a track of particle with mass=pm to a reference plane 
264   // defined by x=xk through media of density=rho and radiationLength=x0
265
266
267   Double_t ymax=AliTOFGeometry::RinTOF()*TMath::Tan(0.5*AliTOFGeometry::GetAlpha());
268   Bool_t skip = kFALSE;
269   Double_t y=GetYat(AliTOFGeometry::RinTOF(),skip);
270   if (skip) {
271     return kFALSE;
272   }
273   if (y > ymax) {
274     if (!Rotate(AliTOFGeometry::GetAlpha())) return kFALSE;
275   } else if (y <-ymax) {
276     if (!Rotate(-AliTOFGeometry::GetAlpha())) return kFALSE;
277   }
278   
279   Double_t x = GetX();
280   Int_t nsteps=Int_t((AliTOFGeometry::Rmin()-x)/0.5); // 0.5 cm Steps
281   for (Int_t istep=0;istep<nsteps;istep++){
282     Float_t xp = x+istep*0.5; 
283     //    GetPropagationParameters(param);  
284     if (!(PropagateTo(xp,0.,0.))) { /* passing 0.,0. as arguments since now
285                                        this method queries TGeo for material budget
286                                     */
287       return kFALSE;
288     }
289     
290   }
291   
292   if (!(PropagateTo(AliTOFGeometry::RinTOF()))) return kFALSE;
293
294   return kTRUE;
295   
296 }     
297
298 //_________________________________________________________________________
299 Double_t AliTOFtrack::GetPredictedChi2(const AliCluster3D *c) const {
300   //
301   // Estimate the chi2 of the space point "c" with its cov. matrix elements
302   //
303
304   Double_t p[3]={c->GetX(), c->GetY(), c->GetZ()};
305   Double_t covyz[3]={c->GetSigmaY2(), c->GetSigmaYZ(), c->GetSigmaZ2()};
306   Double_t covxyz[3]={c->GetSigmaX2(), c->GetSigmaXY(), c->GetSigmaXZ()};
307   return AliExternalTrackParam::GetPredictedChi2(p, covyz, covxyz);
308 }
309 //_________________________________________________________________________
310 Bool_t AliTOFtrack::PropagateTo(const AliCluster3D *c) {
311   //
312   // Propagates a track to the plane containing cluster 'c'
313   //
314   Double_t oldX=GetX(), oldY=GetY(), oldZ=GetZ();
315   Double_t p[3]={c->GetX(), c->GetY(), c->GetZ()};
316   Double_t covyz[3]={c->GetSigmaY2(), c->GetSigmaYZ(), c->GetSigmaZ2()};
317   Double_t covxyz[3]={c->GetSigmaX2(), c->GetSigmaXY(), c->GetSigmaXZ()};
318   Double_t bz=GetBz();
319   if (!AliExternalTrackParam::PropagateTo(p, covyz, covxyz, bz)) return kFALSE;
320   if (IsStartedTimeIntegral()) {
321     Double_t d = TMath::Sqrt((GetX()-oldX)*(GetX()-oldX) +
322                              (GetY()-oldY)*(GetY()-oldY) +
323                              (GetZ()-oldZ)*(GetZ()-oldZ));
324     if (GetX()<oldX) d=-d;
325     AddTimeStep(d);
326
327   }
328
329   if (GetY() >  GetX()*TMath::Tan(0.5*AliTOFGeometry::GetAlpha())) {
330     if (!(Rotate(AliTOFGeometry::GetAlpha()))) return kFALSE;
331   } else if (GetY() < -GetX()*TMath::Tan(0.5*AliTOFGeometry::GetAlpha())) {
332     if (!(Rotate(-AliTOFGeometry::GetAlpha()))) return kFALSE;
333   }
334
335   return kTRUE;
336 }
337 //_________________________________________________________________________
338 Double_t AliTOFtrack::GetYat(Double_t xk, Bool_t & skip) const {     
339 //-----------------------------------------------------------------
340 // This function calculates the Y-coordinate of a track at the plane x=xk.
341 // Needed for matching with the TOF (I.Belikov)
342 //-----------------------------------------------------------------
343      Double_t y=0.;
344      skip=(!GetYAt(xk,GetBz(),y));
345      return y;
346 }
347
348 //_____________________________________________________________________________
349 Int_t AliTOFtrack::Compare(const TObject *o) const {
350   //-----------------------------------------------------------------
351   // This function compares tracks according to the their curvature
352   //-----------------------------------------------------------------
353   AliTOFtrack *t=(AliTOFtrack*)o;
354   Double_t co=t->GetSigmaY2()*t->GetSigmaZ2();
355   Double_t c =GetSigmaY2()*GetSigmaZ2();
356   if (c>co) return 1;
357   else if (c<co) return -1;
358   return 0;
359 }
360
361
362 //_____________________________________________________________________________
363 void AliTOFtrack::GetPropagationParameters(Double_t *param) {
364
365  //Get average medium density, x0 while propagating the track
366
367   //For TRD holes description
368   /*
369   Double_t thetamin = (90.-31.1) * TMath::Pi()/180.;
370   Double_t thetamax = (90.+31.1) * TMath::Pi()/180.;
371
372   Double_t zmin = -55.;
373   Double_t zmax =  55.;
374   */
375
376   // Detector inner/outer radii
377   Double_t rTPC    = 261.53;
378   Double_t rTPCTRD = 294.5;
379   Double_t rTRD    = 369.1;
380
381   // Medium parameters
382   Double_t x0TPC = 40.;
383   Double_t rhoTPC =0.06124;
384
385   Double_t x0Air = 36.66;
386   Double_t rhoAir =1.2931e-3;
387
388   Double_t x0TRD = 171.7;
389   Double_t rhoTRD =0.33;
390
391   //  Int_t isec = GetSector();
392   Double_t r[3]; GetXYZ(r);
393   //  Float_t thetatr = TMath::ATan2(TMath::Sqrt(r[0]*r[0]+r[1]*r[1]),r[2]);
394
395   /*
396   if(holes){
397     if (isec == 0 || isec == 1 || isec == 2 ) {
398       if( thetatr>=thetamin && thetatr<=thetamax){ 
399         x0TRD= x0Air;
400         rhoTRD = rhoAir;
401       }
402     }
403     if (isec == 11 || isec == 12 || isec == 13 || isec == 14 || isec == 15 ) {
404       if( r[2]>=zmin && r[2]<=zmax){ 
405         x0TRD= x0Air;
406         rhoTRD = rhoAir;
407       }
408     }
409   }
410   */
411   if(GetX() <= rTPC)
412     {param[0]=x0TPC;param[1]=rhoTPC;}
413   else if(GetX() > rTPC &&  GetX() < rTPCTRD)
414     {param[0]=x0Air;param[1]=rhoAir;}
415   else if(GetX() >= rTPCTRD &&  GetX() < rTRD)
416     {param[0]=x0TRD;param[1]=rhoTRD;}
417   else if(GetX() >= rTRD )
418     {param[0]=x0Air;param[1]=rhoAir;}
419 }