]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ITS/AliITStrackV2.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Adding MUON HLT code to the repository.
[u/mrichter/AliRoot.git] / ITS / AliITStrackV2.cxx
CommitLineData
006b5f7f 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 ITS track class
18//
19// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch
a9a2d814 20// dEdx analysis by: Boris Batyunya, JINR, Boris.Batiounia@cern.ch
006b5f7f 21//-------------------------------------------------------------------------
006b5f7f 22#include <TMath.h>
006b5f7f 23
24#include "AliCluster.h"
83d73500 25#include "AliESDtrack.h"
006b5f7f 26#include "AliITStrackV2.h"
27
1aedd96e 28ClassImp(AliITStrackV2)
8676d691 29
e43c066c 30const Int_t kWARN=5;
31
22b13da0 32//____________________________________________________________________________
e6ba4672 33AliITStrackV2::AliITStrackV2():AliKalmanTrack(),
22b13da0 34 fX(0),
35 fAlpha(0),
36 fdEdx(0),
37 fP0(0),
38 fP1(0),
39 fP2(0),
40 fP3(0),
41 fP4(0),
42 fC00(0),
43 fC10(0),
44 fC11(0),
45 fC20(0),
46 fC21(0),
47 fC22(0),
48 fC30(0),
49 fC31(0),
50 fC32(0),
51 fC33(0),
52 fC40(0),
53 fC41(0),
54 fC42(0),
55 fC43(0),
83d73500 56 fC44(0),
57 fESDtrack(0)
e43c066c 58{
15dd636f 59 for(Int_t i=0; i<kMaxLayer; i++) fIndex[i]=0;
e43c066c 60 for(Int_t i=0; i<4; i++) fdEdxSample[i]=0;
e43c066c 61}
62
83d73500 63
83d73500 64//____________________________________________________________________________
67c3dcbe 65AliITStrackV2::AliITStrackV2(AliESDtrack& t,Bool_t c) throw (const Char_t *) :
83d73500 66AliKalmanTrack() {
67 //------------------------------------------------------------------
67c3dcbe 68 // Conversion ESD track -> ITS track.
69 // If c==kTRUE, create the ITS track out of the constrained params.
83d73500 70 //------------------------------------------------------------------
71 SetNumberOfClusters(t.GetITSclusters(fIndex));
72 SetLabel(t.GetLabel());
73 SetMass(t.GetMass());
e43c066c 74 //
75 //
83d73500 76
77 fdEdx=t.GetITSsignal();
78 fAlpha = t.GetAlpha();
79 if (fAlpha < -TMath::Pi()) fAlpha += 2*TMath::Pi();
80 else if (fAlpha >= TMath::Pi()) fAlpha -= 2*TMath::Pi();
81
9b280d80 82 //Conversion of the track parameters
67c3dcbe 83 Double_t x,p[5];
84 if (c) t.GetConstrainedExternalParameters(x,p);
85 else t.GetExternalParameters(x,p);
9b280d80 86 fX=x; x=GetConvConst();
006b5f7f 87 fP0=p[0];
88 fP1=p[1];
89 fP2=p[2];
90 fP3=p[3];
91 fP4=p[4]/x;
92
9b280d80 93 //Conversion of the covariance matrix
67c3dcbe 94 Double_t cv[15];
95 if (c) t.GetConstrainedExternalCovariance(cv);
96 else t.GetExternalCovariance(cv);
97 fC00=cv[0 ];
98 fC10=cv[1 ]; fC11=cv[2 ];
99 fC20=cv[3 ]; fC21=cv[4 ]; fC22=cv[5 ];
100 fC30=cv[6 ]; fC31=cv[7 ]; fC32=cv[8 ]; fC33=cv[9 ];
101 fC40=cv[10]/x; fC41=cv[11]/x; fC42=cv[12]/x; fC43=cv[13]/x; fC44=cv[14]/x/x;
006b5f7f 102
83d73500 103 if (t.GetStatus()&AliESDtrack::kTIME) {
104 StartTimeIntegral();
105 Double_t times[10]; t.GetIntegratedTimes(times); SetIntegratedTimes(times);
106 SetIntegratedLength(t.GetIntegratedLength());
107 }
108 fESDtrack=&t;
e43c066c 109
81e97e0d 110 // if (!Invariant()) throw "AliITStrackV2: conversion failed !\n";
68b8060b 111 for(Int_t i=0; i<4; i++) fdEdxSample[i]=0;
112
006b5f7f 113}
114
15dd636f 115void AliITStrackV2::UpdateESDtrack(ULong_t flags) const {
83d73500 116 fESDtrack->UpdateTrackParams(this,flags);
117}
15dd636f 118void AliITStrackV2::SetConstrainedESDtrack(Double_t chi2) const {
67c3dcbe 119 fESDtrack->SetConstrainedTrackParams(this,chi2);
120}
83d73500 121
006b5f7f 122//____________________________________________________________________________
123AliITStrackV2::AliITStrackV2(const AliITStrackV2& t) : AliKalmanTrack(t) {
124 //------------------------------------------------------------------
125 //Copy constructor
126 //------------------------------------------------------------------
127 fX=t.fX;
128 fAlpha=t.fAlpha;
129 fdEdx=t.fdEdx;
130
131 fP0=t.fP0; fP1=t.fP1; fP2=t.fP2; fP3=t.fP3; fP4=t.fP4;
132
133 fC00=t.fC00;
134 fC10=t.fC10; fC11=t.fC11;
135 fC20=t.fC20; fC21=t.fC21; fC22=t.fC22;
136 fC30=t.fC30; fC31=t.fC31; fC32=t.fC32; fC33=t.fC33;
137 fC40=t.fC40; fC41=t.fC41; fC42=t.fC42; fC43=t.fC43; fC44=t.fC44;
138
139 Int_t n=GetNumberOfClusters();
23efe5f1 140 for (Int_t i=0; i<n; i++) {
e43c066c 141 fIndex[i]=t.fIndex[i];
142 if (i<4) fdEdxSample[i]=t.fdEdxSample[i];
23efe5f1 143 }
83d73500 144 fESDtrack=t.fESDtrack;
006b5f7f 145}
a9a2d814 146
006b5f7f 147//_____________________________________________________________________________
148Int_t AliITStrackV2::Compare(const TObject *o) const {
149 //-----------------------------------------------------------------
150 // This function compares tracks according to the their curvature
151 //-----------------------------------------------------------------
7f6ddf58 152 AliITStrackV2 *t=(AliITStrackV2*)o;
a9a2d814 153 //Double_t co=TMath::Abs(t->Get1Pt());
154 //Double_t c =TMath::Abs(Get1Pt());
15dd636f 155 Double_t co=t->GetSigmaY2()*t->GetSigmaZ2();
156 Double_t c =GetSigmaY2()*GetSigmaZ2();
006b5f7f 157 if (c>co) return 1;
158 else if (c<co) return -1;
159 return 0;
160}
161
162//_____________________________________________________________________________
163void AliITStrackV2::GetExternalCovariance(Double_t cc[15]) const {
164 //-------------------------------------------------------------------------
165 // This function returns an external representation of the covriance matrix.
166 // (See comments in AliTPCtrack.h about external track representation)
167 //-------------------------------------------------------------------------
9b280d80 168 Double_t a=GetConvConst();
006b5f7f 169
170 cc[0 ]=fC00;
171 cc[1 ]=fC10; cc[2 ]=fC11;
172 cc[3 ]=fC20; cc[4 ]=fC21; cc[5 ]=fC22;
173 cc[6 ]=fC30; cc[7 ]=fC31; cc[8 ]=fC32; cc[9 ]=fC33;
174 cc[10]=fC40*a; cc[11]=fC41*a; cc[12]=fC42*a; cc[13]=fC43*a; cc[14]=fC44*a*a;
175}
176
177//____________________________________________________________________________
a9a2d814 178Int_t AliITStrackV2::PropagateToVertex(Double_t d,Double_t x0) {
006b5f7f 179 //------------------------------------------------------------------
180 //This function propagates a track to the minimal distance from the origin
181 //------------------------------------------------------------------
5b316485 182 //Double_t xv=fP2*(fX*fP2 - fP0*TMath::Sqrt(1.- fP2*fP2)); //linear approxim.
183 Double_t tgf=-(fP4*fX - fP2)/(fP4*fP0 + TMath::Sqrt(1 - fP2*fP2));
184 Double_t snf=tgf/TMath::Sqrt(1.+ tgf*tgf);
185 Double_t xv=(snf - fP2)/fP4 + fX;
186 return PropagateTo(xv,d,x0);
006b5f7f 187}
188
189//____________________________________________________________________________
190Int_t AliITStrackV2::
191GetGlobalXYZat(Double_t xk, Double_t &x, Double_t &y, Double_t &z) const {
192 //------------------------------------------------------------------
193 //This function returns a track position in the global system
194 //------------------------------------------------------------------
195 Double_t dx=xk-fX;
196 Double_t f1=fP2, f2=f1 + fP4*dx;
a9a2d814 197 if (TMath::Abs(f2) >= 0.9999) {
7f6ddf58 198 Int_t n=GetNumberOfClusters();
199 if (n>kWARN)
8676d691 200 Warning("GetGlobalXYZat","Propagation failed (%d) !\n",n);
006b5f7f 201 return 0;
202 }
203
204 Double_t r1=sqrt(1.- f1*f1), r2=sqrt(1.- f2*f2);
205
206 Double_t yk = fP0 + dx*(f1+f2)/(r1+r2);
207 Double_t zk = fP1 + dx*(f1+f2)/(f1*r2 + f2*r1)*fP3;
208
209 Double_t cs=TMath::Cos(fAlpha), sn=TMath::Sin(fAlpha);
210 x = xk*cs - yk*sn;
211 y = xk*sn + yk*cs;
212 z = zk;
213
214 return 1;
215}
216
217//_____________________________________________________________________________
218Double_t AliITStrackV2::GetPredictedChi2(const AliCluster *c) const
219{
220 //-----------------------------------------------------------------
221 // This function calculates a predicted chi2 increment.
222 //-----------------------------------------------------------------
223 Double_t r00=c->GetSigmaY2(), r01=0., r11=c->GetSigmaZ2();
224 r00+=fC00; r01+=fC10; r11+=fC11;
babd135a 225 //
006b5f7f 226 Double_t det=r00*r11 - r01*r01;
7f6ddf58 227 if (TMath::Abs(det) < 1.e-30) {
006b5f7f 228 Int_t n=GetNumberOfClusters();
7f6ddf58 229 if (n>kWARN)
8676d691 230 Warning("GetPredictedChi2","Singular matrix (%d) !\n",n);
006b5f7f 231 return 1e10;
232 }
233 Double_t tmp=r00; r00=r11; r11=tmp; r01=-r01;
33329428 234
006b5f7f 235 Double_t dy=c->GetY() - fP0, dz=c->GetZ() - fP1;
33329428 236
006b5f7f 237 return (dy*r00*dy + 2*r01*dy*dz + dz*r11*dz)/det;
238}
239
006b5f7f 240//____________________________________________________________________________
a9a2d814 241Int_t AliITStrackV2::CorrectForMaterial(Double_t d, Double_t x0) {
242 //------------------------------------------------------------------
243 //This function corrects the track parameters for crossed material
244 //------------------------------------------------------------------
e43c066c 245 Double_t p2=(1.+ fP3*fP3)/(Get1Pt()*Get1Pt());
15dd636f 246 Double_t beta2=p2/(p2 + GetMass()*GetMass());
a9a2d814 247 d*=TMath::Sqrt((1.+ fP3*fP3)/(1.- fP2*fP2));
248
249 //Multiple scattering******************
250 if (d!=0) {
33329428 251 Double_t theta2=14.1*14.1/(beta2*p2*1e6)*TMath::Abs(d);
252 //Double_t theta2=1.0259e-6*14*14/28/(beta2*p2)*TMath::Abs(d)*9.36*2.33;
a9a2d814 253 fC22 += theta2*(1.- fP2*fP2)*(1. + fP3*fP3);
254 fC33 += theta2*(1. + fP3*fP3)*(1. + fP3*fP3);
255 fC43 += theta2*fP3*fP4*(1. + fP3*fP3);
256 fC44 += theta2*fP3*fP4*fP3*fP4;
257 }
258
259 //Energy losses************************
260 if (x0!=0.) {
261 d*=x0;
15dd636f 262 Double_t dE=0.153e-3/beta2*(log(5940*beta2/(1-beta2)) - beta2)*d;
263 if (beta2/(1-beta2)>3.5*3.5)
83d73500 264 dE=0.153e-3/beta2*(log(3.5*5940)+0.5*log(beta2/(1-beta2)) - beta2)*d;
15dd636f 265
266 fP4*=(1.- TMath::Sqrt(p2+GetMass()*GetMass())/p2*dE);
a9a2d814 267 }
268
269 if (!Invariant()) return 0;
270
271 return 1;
272}
273
274//____________________________________________________________________________
275Int_t AliITStrackV2::PropagateTo(Double_t xk, Double_t d, Double_t x0) {
006b5f7f 276 //------------------------------------------------------------------
277 //This function propagates a track
278 //------------------------------------------------------------------
279 Double_t x1=fX, x2=xk, dx=x2-x1;
280 Double_t f1=fP2, f2=f1 + fP4*dx;
81e97e0d 281 if (TMath::Abs(f2) >= 0.98) {
282 // MI change - don't propagate highly inclined tracks
283 // covariance matrix distorted
284 // Int_t n=GetNumberOfClusters();
285 // if (n>kWARN)
286 // Warning("PropagateTo","Propagation failed !\n",n);
006b5f7f 287 return 0;
288 }
289
f29dbb4b 290 // old position [SR, GSI, 17.02.2003]
291 Double_t oldX = fX, oldY = fP0, oldZ = fP1;
292
006b5f7f 293 Double_t r1=sqrt(1.- f1*f1), r2=sqrt(1.- f2*f2);
294
295 fP0 += dx*(f1+f2)/(r1+r2);
296 fP1 += dx*(f1+f2)/(f1*r2 + f2*r1)*fP3;
297 fP2 += dx*fP4;
298
299 //f = F - 1
300
301 Double_t f02= dx/(r1*r1*r1);
302 Double_t f04=0.5*dx*dx/(r1*r1*r1);
303 Double_t f12= dx*fP3*f1/(r1*r1*r1);
304 Double_t f14=0.5*dx*dx*fP3*f1/(r1*r1*r1);
305 Double_t f13= dx/r1;
306 Double_t f24= dx;
307
308 //b = C*ft
309 Double_t b00=f02*fC20 + f04*fC40, b01=f12*fC20 + f14*fC40 + f13*fC30;
310 Double_t b02=f24*fC40;
311 Double_t b10=f02*fC21 + f04*fC41, b11=f12*fC21 + f14*fC41 + f13*fC31;
312 Double_t b12=f24*fC41;
313 Double_t b20=f02*fC22 + f04*fC42, b21=f12*fC22 + f14*fC42 + f13*fC32;
314 Double_t b22=f24*fC42;
315 Double_t b40=f02*fC42 + f04*fC44, b41=f12*fC42 + f14*fC44 + f13*fC43;
316 Double_t b42=f24*fC44;
317 Double_t b30=f02*fC32 + f04*fC43, b31=f12*fC32 + f14*fC43 + f13*fC33;
318 Double_t b32=f24*fC43;
319
320 //a = f*b = f*C*ft
321 Double_t a00=f02*b20+f04*b40,a01=f02*b21+f04*b41,a02=f02*b22+f04*b42;
322 Double_t a11=f12*b21+f14*b41+f13*b31,a12=f12*b22+f14*b42+f13*b32;
323 Double_t a22=f24*b42;
324
325 //F*C*Ft = C + (b + bt + a)
326 fC00 += b00 + b00 + a00;
327 fC10 += b10 + b01 + a01;
328 fC20 += b20 + b02 + a02;
329 fC30 += b30;
330 fC40 += b40;
331 fC11 += b11 + b11 + a11;
332 fC21 += b21 + b12 + a12;
333 fC31 += b31;
334 fC41 += b41;
335 fC22 += b22 + b22 + a22;
336 fC32 += b32;
337 fC42 += b42;
338
339 fX=x2;
340
a9a2d814 341 if (!CorrectForMaterial(d,x0)) return 0;
006b5f7f 342
f29dbb4b 343 // Integrated Time [SR, GSI, 17.02.2003]
880f41b9 344 if (IsStartedTimeIntegral() && fX>oldX) {
8676d691 345 Double_t l2 = (fX-oldX)*(fX-oldX)+(fP0-oldY)*(fP0-oldY)+
346 (fP1-oldZ)*(fP1-oldZ);
f29dbb4b 347 AddTimeStep(TMath::Sqrt(l2));
348 }
349 //
350
006b5f7f 351 return 1;
352}
353
354//____________________________________________________________________________
355Int_t AliITStrackV2::Update(const AliCluster* c, Double_t chi2, UInt_t index) {
356 //------------------------------------------------------------------
357 //This function updates track parameters
358 //------------------------------------------------------------------
359 Double_t p0=fP0,p1=fP1,p2=fP2,p3=fP3,p4=fP4;
360 Double_t c00=fC00;
361 Double_t c10=fC10, c11=fC11;
362 Double_t c20=fC20, c21=fC21, c22=fC22;
363 Double_t c30=fC30, c31=fC31, c32=fC32, c33=fC33;
364 Double_t c40=fC40, c41=fC41, c42=fC42, c43=fC43, c44=fC44;
365
366
367 Double_t r00=c->GetSigmaY2(), r01=0., r11=c->GetSigmaZ2();
368 r00+=fC00; r01+=fC10; r11+=fC11;
369 Double_t det=r00*r11 - r01*r01;
370 Double_t tmp=r00; r00=r11/det; r11=tmp/det; r01=-r01/det;
371
babd135a 372
006b5f7f 373 Double_t k00=fC00*r00+fC10*r01, k01=fC00*r01+fC10*r11;
374 Double_t k10=fC10*r00+fC11*r01, k11=fC10*r01+fC11*r11;
375 Double_t k20=fC20*r00+fC21*r01, k21=fC20*r01+fC21*r11;
376 Double_t k30=fC30*r00+fC31*r01, k31=fC30*r01+fC31*r11;
377 Double_t k40=fC40*r00+fC41*r01, k41=fC40*r01+fC41*r11;
378
379 Double_t dy=c->GetY() - fP0, dz=c->GetZ() - fP1;
380 Double_t sf=fP2 + k20*dy + k21*dz;
7f6ddf58 381
006b5f7f 382 fP0 += k00*dy + k01*dz;
383 fP1 += k10*dy + k11*dz;
384 fP2 = sf;
385 fP3 += k30*dy + k31*dz;
386 fP4 += k40*dy + k41*dz;
387
388 Double_t c01=fC10, c02=fC20, c03=fC30, c04=fC40;
389 Double_t c12=fC21, c13=fC31, c14=fC41;
390
391 fC00-=k00*fC00+k01*fC10; fC10-=k00*c01+k01*fC11;
392 fC20-=k00*c02+k01*c12; fC30-=k00*c03+k01*c13;
393 fC40-=k00*c04+k01*c14;
394
395 fC11-=k10*c01+k11*fC11;
396 fC21-=k10*c02+k11*c12; fC31-=k10*c03+k11*c13;
397 fC41-=k10*c04+k11*c14;
398
399 fC22-=k20*c02+k21*c12; fC32-=k20*c03+k21*c13;
400 fC42-=k20*c04+k21*c14;
401
402 fC33-=k30*c03+k31*c13;
403 fC43-=k30*c04+k31*c14;
404
405 fC44-=k40*c04+k41*c14;
406
407 if (!Invariant()) {
408 fP0=p0; fP1=p1; fP2=p2; fP3=p3; fP4=p4;
409 fC00=c00;
410 fC10=c10; fC11=c11;
411 fC20=c20; fC21=c21; fC22=c22;
412 fC30=c30; fC31=c31; fC32=c32; fC33=c33;
413 fC40=c40; fC41=c41; fC42=c42; fC43=c43; fC44=c44;
414 return 0;
415 }
416
67c3dcbe 417 if (chi2<0) return 1;
418
006b5f7f 419 Int_t n=GetNumberOfClusters();
420 fIndex[n]=index;
421 SetNumberOfClusters(n+1);
422 SetChi2(GetChi2()+chi2);
423
424 return 1;
425}
426
006b5f7f 427Int_t AliITStrackV2::Invariant() const {
428 //------------------------------------------------------------------
429 // This function is for debugging purpose only
430 //------------------------------------------------------------------
7f6ddf58 431 Int_t n=GetNumberOfClusters();
432
a9a2d814 433 if (TMath::Abs(fP2)>=0.9999){
8676d691 434 if (n>kWARN) Warning("Invariant","fP2=%f\n",fP2);
a9a2d814 435 return 0;
436 }
437 if (fC00<=0 || fC00>9.) {
8676d691 438 if (n>kWARN) Warning("Invariant","fC00=%f\n",fC00);
a9a2d814 439 return 0;
440 }
441 if (fC11<=0 || fC11>9.) {
8676d691 442 if (n>kWARN) Warning("Invariant","fC11=%f\n",fC11);
a9a2d814 443 return 0;
444 }
445 if (fC22<=0 || fC22>1.) {
8676d691 446 if (n>kWARN) Warning("Invariant","fC22=%f\n",fC22);
a9a2d814 447 return 0;
448 }
449 if (fC33<=0 || fC33>1.) {
8676d691 450 if (n>kWARN) Warning("Invariant","fC33=%f\n",fC33);
a9a2d814 451 return 0;
452 }
453 if (fC44<=0 || fC44>6e-5) {
8676d691 454 if (n>kWARN) Warning("Invariant","fC44=%f\n",fC44);
a9a2d814 455 return 0;
14825d5a 456 }
14825d5a 457 return 1;
006b5f7f 458}
459
460//____________________________________________________________________________
a9a2d814 461Int_t AliITStrackV2::Propagate(Double_t alp,Double_t xk) {
006b5f7f 462 //------------------------------------------------------------------
463 //This function propagates a track
464 //------------------------------------------------------------------
33329428 465 Double_t alpha=fAlpha, x=fX;
006b5f7f 466 Double_t p0=fP0,p1=fP1,p2=fP2,p3=fP3,p4=fP4;
467 Double_t c00=fC00;
468 Double_t c10=fC10, c11=fC11;
469 Double_t c20=fC20, c21=fC21, c22=fC22;
470 Double_t c30=fC30, c31=fC31, c32=fC32, c33=fC33;
471 Double_t c40=fC40, c41=fC41, c42=fC42, c43=fC43, c44=fC44;
472
33329428 473 if (alp < -TMath::Pi()) alp += 2*TMath::Pi();
474 else if (alp >= TMath::Pi()) alp -= 2*TMath::Pi();
475 Double_t ca=TMath::Cos(alp-fAlpha), sa=TMath::Sin(alp-fAlpha);
476 Double_t sf=fP2, cf=TMath::Sqrt(1.- fP2*fP2);
006b5f7f 477
33329428 478 // **** rotation **********************
479 {
006b5f7f 480 fAlpha = alp;
33329428 481 fX = x*ca + p0*sa;
482 fP0= -x*sa + p0*ca;
483 fP2= sf*ca - cf*sa;
484
c219fd63 485 Double_t rr=(ca+sf/cf*sa);
486
487 fC00 *= (ca*ca);
488 fC10 *= ca;
489 fC20 *= ca*rr;
490 fC30 *= ca;
491 fC40 *= ca;
492 //fC11 = fC11;
493 fC21 *= rr;
494 //fC31 = fC31;
495 //fC41 = fC41;
496 fC22 *= rr*rr;
497 fC32 *= rr;
498 fC42 *= rr;
499 //fC33=fC33;
500 //fC43=fC43;
501 //fC44=fC44;
502
33329428 503 }
006b5f7f 504
33329428 505 // **** translation ******************
506 {
507 Double_t dx=xk-fX;
006b5f7f 508 Double_t f1=fP2, f2=f1 + fP4*dx;
81e97e0d 509 if (TMath::Abs(f2) >= 0.98) {
510 // don't propagate highly inclined tracks MI
006b5f7f 511 return 0;
512 }
81e97e0d 513 // if (TMath::Abs(f2) >= 0.9999) {
514// Int_t n=GetNumberOfClusters();
515// if (n>kWARN)
516// Warning("Propagate","Propagation failed (%d) !\n",n);
517// return 0;
518// }
33329428 519 Double_t r1=TMath::Sqrt(1.- f1*f1), r2=TMath::Sqrt(1.- f2*f2);
006b5f7f 520
33329428 521 fX=xk;
006b5f7f 522 fP0 += dx*(f1+f2)/(r1+r2);
523 fP1 += dx*(f1+f2)/(f1*r2 + f2*r1)*fP3;
524 fP2 += dx*fP4;
525
c219fd63 526 //f = F - 1
527
528 Double_t f02= dx/(r1*r1*r1);
529 Double_t f04=0.5*dx*dx/(r1*r1*r1);
530 Double_t f12= dx*fP3*f1/(r1*r1*r1);
531 Double_t f14=0.5*dx*dx*fP3*f1/(r1*r1*r1);
532 Double_t f13= dx/r1;
533 Double_t f24= dx;
534
535 //b = C*ft
536 Double_t b00=f02*fC20 + f04*fC40, b01=f12*fC20 + f14*fC40 + f13*fC30;
537 Double_t b02=f24*fC40;
538 Double_t b10=f02*fC21 + f04*fC41, b11=f12*fC21 + f14*fC41 + f13*fC31;
539 Double_t b12=f24*fC41;
540 Double_t b20=f02*fC22 + f04*fC42, b21=f12*fC22 + f14*fC42 + f13*fC32;
541 Double_t b22=f24*fC42;
542 Double_t b40=f02*fC42 + f04*fC44, b41=f12*fC42 + f14*fC44 + f13*fC43;
543 Double_t b42=f24*fC44;
544 Double_t b30=f02*fC32 + f04*fC43, b31=f12*fC32 + f14*fC43 + f13*fC33;
545 Double_t b32=f24*fC43;
546
547 //a = f*b = f*C*ft
548 Double_t a00=f02*b20+f04*b40,a01=f02*b21+f04*b41,a02=f02*b22+f04*b42;
549 Double_t a11=f12*b21+f14*b41+f13*b31,a12=f12*b22+f14*b42+f13*b32;
550 Double_t a22=f24*b42;
551
552 //F*C*Ft = C + (b + bt + a)
553 fC00 += b00 + b00 + a00;
554 fC10 += b10 + b01 + a01;
555 fC20 += b20 + b02 + a02;
556 fC30 += b30;
557 fC40 += b40;
558 fC11 += b11 + b11 + a11;
559 fC21 += b21 + b12 + a12;
560 fC31 += b31;
561 fC41 += b41;
562 fC22 += b22 + b22 + a22;
563 fC32 += b32;
564 fC42 += b42;
006b5f7f 565
006b5f7f 566 if (!Invariant()) {
33329428 567 fAlpha=alpha;
568 fX=x;
006b5f7f 569 fP0=p0; fP1=p1; fP2=p2; fP3=p3; fP4=p4;
570 fC00=c00;
571 fC10=c10; fC11=c11;
572 fC20=c20; fC21=c21; fC22=c22;
573 fC30=c30; fC31=c31; fC32=c32; fC33=c33;
574 fC40=c40; fC41=c41; fC42=c42; fC43=c43; fC44=c44;
575 return 0;
576 }
33329428 577 }
006b5f7f 578
579 return 1;
580}
581
c219fd63 582
a6a6f799 583Double_t AliITStrackV2::GetD(Double_t x, Double_t y) const {
006b5f7f 584 //------------------------------------------------------------------
a6a6f799 585 // This function calculates the transverse impact parameter
586 // with respect to a point with global coordinates (x,y)
006b5f7f 587 //------------------------------------------------------------------
a6a6f799 588 Double_t xt=fX, yt=fP0;
589
590 Double_t sn=TMath::Sin(fAlpha), cs=TMath::Cos(fAlpha);
591 Double_t a = x*cs + y*sn;
592 y = -x*sn + y*cs; x=a;
593 xt-=x; yt-=y;
594
595 sn=fP4*xt - fP2; cs=fP4*yt + TMath::Sqrt(1.- fP2*fP2);
596 a=2*(xt*fP2 - yt*TMath::Sqrt(1.- fP2*fP2))-fP4*(xt*xt + yt*yt);
006b5f7f 597 if (fP4<0) a=-a;
598 return a/(1 + TMath::Sqrt(sn*sn + cs*cs));
a9a2d814 599}
006b5f7f 600
babd135a 601Double_t AliITStrackV2::GetZat(Double_t x) const {
602 //------------------------------------------------------------------
603 // This function calculates the z at given x point - in current coordinate system
604 //------------------------------------------------------------------
605 Double_t x1=fX, x2=x, dx=x2-x1;
606 //
607 Double_t f1=fP2, f2=f1 + fP4*dx;
608 if (TMath::Abs(f2) >= 0.9999) {
609 return 10000000;
610 }
babd135a 611 Double_t r1=sqrt(1.- f1*f1), r2=sqrt(1.- f2*f2);
612 Double_t z = fP1 + dx*(f1+f2)/(f1*r2 + f2*r1)*fP3;
613 return z;
614}
615
616
617
618
8676d691 619Int_t AliITStrackV2::Improve(Double_t x0,Double_t xyz[3],Double_t ers[3]) {
a9a2d814 620 //------------------------------------------------------------------
621 //This function improves angular track parameters
622 //------------------------------------------------------------------
8676d691 623 Double_t cs=TMath::Cos(fAlpha), sn=TMath::Sin(fAlpha);
624 //Double_t xv = xyz[0]*cs + xyz[1]*sn; // vertex
625 Double_t yv =-xyz[0]*sn + xyz[1]*cs; // in the
626 Double_t zv = xyz[2]; // local frame
a9a2d814 627 Double_t dy=fP0-yv, dz=fP1-zv;
628 Double_t r2=fX*fX+dy*dy;
629 Double_t p2=(1.+ GetTgl()*GetTgl())/(Get1Pt()*Get1Pt());
630 Double_t beta2=p2/(p2 + GetMass()*GetMass());
631 x0*=TMath::Sqrt((1.+ GetTgl()*GetTgl())/(1.- GetSnp()*GetSnp()));
8676d691 632 Double_t theta2=14.1*14.1/(beta2*p2*1e6)*x0;
633 //Double_t theta2=1.0259e-6*14*14/28/(beta2*p2)*x0*9.36*2.33;
a9a2d814 634 {
25d954aa 635 Double_t dummy=4/r2-fP4*fP4;
636 if (dummy < 0) return 0;
637 Double_t parp=0.5*(fP4*fX + dy*TMath::Sqrt(dummy));
a9a2d814 638 Double_t sigma2p = theta2*(1.- GetSnp()*GetSnp())*(1. + GetTgl()*GetTgl());
639 sigma2p += fC00/r2*(1.- dy*dy/r2)*(1.- dy*dy/r2);
8676d691 640 sigma2p += ers[1]*ers[1]/r2;
a9a2d814 641 sigma2p += 0.25*fC44*fX*fX;
642 Double_t eps2p=sigma2p/(fC22+sigma2p);
643 fP0 += fC20/(fC22+sigma2p)*(parp-fP2);
644 fP2 = eps2p*fP2 + (1-eps2p)*parp;
645 fC22 *= eps2p;
646 fC20 *= eps2p;
647 }
648 {
649 Double_t parl=0.5*fP4*dz/TMath::ASin(0.5*fP4*TMath::Sqrt(r2));
650 Double_t sigma2l=theta2;
651 sigma2l += fC11/r2+fC00*dy*dy*dz*dz/(r2*r2*r2);
8676d691 652 sigma2l += ers[2]*ers[2]/r2;
a9a2d814 653 Double_t eps2l=sigma2l/(fC33+sigma2l);
654 fP1 += fC31/(fC33+sigma2l)*(parl-fP3);
655 fP4 += fC43/(fC33+sigma2l)*(parl-fP3);
656 fP3 = eps2l*fP3 + (1-eps2l)*parl;
657 fC33 *= eps2l; fC43 *= eps2l;
658 fC31 *= eps2l;
659 }
660 if (!Invariant()) return 0;
661 return 1;
662}
006b5f7f 663
14825d5a 664void AliITStrackV2::ResetCovariance() {
665 //------------------------------------------------------------------
666 //This function makes a track forget its history :)
667 //------------------------------------------------------------------
668
669 fC00*=10.;
670 fC10=0.; fC11*=10.;
671 fC20=0.; fC21=0.; fC22*=10.;
672 fC30=0.; fC31=0.; fC32=0.; fC33*=10.;
673 fC40=0.; fC41=0.; fC42=0.; fC43=0.; fC44*=10.;
674
675}
23efe5f1 676
677void AliITStrackV2::CookdEdx(Double_t low, Double_t up) {
678 //-----------------------------------------------------------------
a9a2d814 679 // This function calculates dE/dX within the "low" and "up" cuts.
680 // Origin: Boris Batyunya, JINR, Boris.Batiounia@cern.ch
23efe5f1 681 //-----------------------------------------------------------------
23efe5f1 682 // The clusters order is: SSD-2, SSD-1, SDD-2, SDD-1, SPD-2, SPD-1
c0dd5278 683
684 Int_t i;
685 Int_t nc=0;
686 for (i=0; i<GetNumberOfClusters(); i++) {
687 Int_t idx=GetClusterIndex(i);
688 idx=(idx&0xf0000000)>>28;
689 if (idx>1) nc++; // Take only SSD and SDD
690 }
23efe5f1 691
a9a2d814 692 Int_t swap;//stupid sorting
23efe5f1 693 do {
694 swap=0;
695 for (i=0; i<nc-1; i++) {
696 if (fdEdxSample[i]<=fdEdxSample[i+1]) continue;
697 Float_t tmp=fdEdxSample[i];
698 fdEdxSample[i]=fdEdxSample[i+1]; fdEdxSample[i+1]=tmp;
699 swap++;
700 }
701 } while (swap);
702
a9a2d814 703 Int_t nl=Int_t(low*nc), nu=Int_t(up*nc); //b.b. to take two lowest dEdX
704 // values from four ones choose
705 // nu=2
23efe5f1 706 Float_t dedx=0;
a9a2d814 707 for (i=nl; i<nu; i++) dedx += fdEdxSample[i];
c0dd5278 708 if (nu-nl>0) dedx /= (nu-nl);
23efe5f1 709
23efe5f1 710 SetdEdx(dedx);
711}