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