]> 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
Initialize arrays in constructor (Ivana)
[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
68b8060b 110 for(Int_t i=0; i<4; i++) fdEdxSample[i]=0;
111
15dd636f 112 if (!Invariant()) throw "AliITStrackV2: conversion failed !\n";
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;
a9a2d814 281 if (TMath::Abs(f2) >= 0.9999) {
006b5f7f 282 Int_t n=GetNumberOfClusters();
7f6ddf58 283 if (n>kWARN)
8676d691 284 Warning("PropagateTo","Propagation failed !\n",n);
006b5f7f 285 return 0;
286 }
287
f29dbb4b 288 // old position [SR, GSI, 17.02.2003]
289 Double_t oldX = fX, oldY = fP0, oldZ = fP1;
290
006b5f7f 291 Double_t r1=sqrt(1.- f1*f1), r2=sqrt(1.- f2*f2);
292
293 fP0 += dx*(f1+f2)/(r1+r2);
294 fP1 += dx*(f1+f2)/(f1*r2 + f2*r1)*fP3;
295 fP2 += dx*fP4;
296
297 //f = F - 1
298
299 Double_t f02= dx/(r1*r1*r1);
300 Double_t f04=0.5*dx*dx/(r1*r1*r1);
301 Double_t f12= dx*fP3*f1/(r1*r1*r1);
302 Double_t f14=0.5*dx*dx*fP3*f1/(r1*r1*r1);
303 Double_t f13= dx/r1;
304 Double_t f24= dx;
305
306 //b = C*ft
307 Double_t b00=f02*fC20 + f04*fC40, b01=f12*fC20 + f14*fC40 + f13*fC30;
308 Double_t b02=f24*fC40;
309 Double_t b10=f02*fC21 + f04*fC41, b11=f12*fC21 + f14*fC41 + f13*fC31;
310 Double_t b12=f24*fC41;
311 Double_t b20=f02*fC22 + f04*fC42, b21=f12*fC22 + f14*fC42 + f13*fC32;
312 Double_t b22=f24*fC42;
313 Double_t b40=f02*fC42 + f04*fC44, b41=f12*fC42 + f14*fC44 + f13*fC43;
314 Double_t b42=f24*fC44;
315 Double_t b30=f02*fC32 + f04*fC43, b31=f12*fC32 + f14*fC43 + f13*fC33;
316 Double_t b32=f24*fC43;
317
318 //a = f*b = f*C*ft
319 Double_t a00=f02*b20+f04*b40,a01=f02*b21+f04*b41,a02=f02*b22+f04*b42;
320 Double_t a11=f12*b21+f14*b41+f13*b31,a12=f12*b22+f14*b42+f13*b32;
321 Double_t a22=f24*b42;
322
323 //F*C*Ft = C + (b + bt + a)
324 fC00 += b00 + b00 + a00;
325 fC10 += b10 + b01 + a01;
326 fC20 += b20 + b02 + a02;
327 fC30 += b30;
328 fC40 += b40;
329 fC11 += b11 + b11 + a11;
330 fC21 += b21 + b12 + a12;
331 fC31 += b31;
332 fC41 += b41;
333 fC22 += b22 + b22 + a22;
334 fC32 += b32;
335 fC42 += b42;
336
337 fX=x2;
338
a9a2d814 339 if (!CorrectForMaterial(d,x0)) return 0;
006b5f7f 340
f29dbb4b 341 // Integrated Time [SR, GSI, 17.02.2003]
880f41b9 342 if (IsStartedTimeIntegral() && fX>oldX) {
8676d691 343 Double_t l2 = (fX-oldX)*(fX-oldX)+(fP0-oldY)*(fP0-oldY)+
344 (fP1-oldZ)*(fP1-oldZ);
f29dbb4b 345 AddTimeStep(TMath::Sqrt(l2));
346 }
347 //
348
006b5f7f 349 return 1;
350}
351
352//____________________________________________________________________________
353Int_t AliITStrackV2::Update(const AliCluster* c, Double_t chi2, UInt_t index) {
354 //------------------------------------------------------------------
355 //This function updates track parameters
356 //------------------------------------------------------------------
357 Double_t p0=fP0,p1=fP1,p2=fP2,p3=fP3,p4=fP4;
358 Double_t c00=fC00;
359 Double_t c10=fC10, c11=fC11;
360 Double_t c20=fC20, c21=fC21, c22=fC22;
361 Double_t c30=fC30, c31=fC31, c32=fC32, c33=fC33;
362 Double_t c40=fC40, c41=fC41, c42=fC42, c43=fC43, c44=fC44;
363
364
365 Double_t r00=c->GetSigmaY2(), r01=0., r11=c->GetSigmaZ2();
366 r00+=fC00; r01+=fC10; r11+=fC11;
367 Double_t det=r00*r11 - r01*r01;
368 Double_t tmp=r00; r00=r11/det; r11=tmp/det; r01=-r01/det;
369
babd135a 370
006b5f7f 371 Double_t k00=fC00*r00+fC10*r01, k01=fC00*r01+fC10*r11;
372 Double_t k10=fC10*r00+fC11*r01, k11=fC10*r01+fC11*r11;
373 Double_t k20=fC20*r00+fC21*r01, k21=fC20*r01+fC21*r11;
374 Double_t k30=fC30*r00+fC31*r01, k31=fC30*r01+fC31*r11;
375 Double_t k40=fC40*r00+fC41*r01, k41=fC40*r01+fC41*r11;
376
377 Double_t dy=c->GetY() - fP0, dz=c->GetZ() - fP1;
378 Double_t sf=fP2 + k20*dy + k21*dz;
7f6ddf58 379
006b5f7f 380 fP0 += k00*dy + k01*dz;
381 fP1 += k10*dy + k11*dz;
382 fP2 = sf;
383 fP3 += k30*dy + k31*dz;
384 fP4 += k40*dy + k41*dz;
385
386 Double_t c01=fC10, c02=fC20, c03=fC30, c04=fC40;
387 Double_t c12=fC21, c13=fC31, c14=fC41;
388
389 fC00-=k00*fC00+k01*fC10; fC10-=k00*c01+k01*fC11;
390 fC20-=k00*c02+k01*c12; fC30-=k00*c03+k01*c13;
391 fC40-=k00*c04+k01*c14;
392
393 fC11-=k10*c01+k11*fC11;
394 fC21-=k10*c02+k11*c12; fC31-=k10*c03+k11*c13;
395 fC41-=k10*c04+k11*c14;
396
397 fC22-=k20*c02+k21*c12; fC32-=k20*c03+k21*c13;
398 fC42-=k20*c04+k21*c14;
399
400 fC33-=k30*c03+k31*c13;
401 fC43-=k30*c04+k31*c14;
402
403 fC44-=k40*c04+k41*c14;
404
405 if (!Invariant()) {
406 fP0=p0; fP1=p1; fP2=p2; fP3=p3; fP4=p4;
407 fC00=c00;
408 fC10=c10; fC11=c11;
409 fC20=c20; fC21=c21; fC22=c22;
410 fC30=c30; fC31=c31; fC32=c32; fC33=c33;
411 fC40=c40; fC41=c41; fC42=c42; fC43=c43; fC44=c44;
412 return 0;
413 }
414
67c3dcbe 415 if (chi2<0) return 1;
416
006b5f7f 417 Int_t n=GetNumberOfClusters();
418 fIndex[n]=index;
419 SetNumberOfClusters(n+1);
420 SetChi2(GetChi2()+chi2);
421
422 return 1;
423}
424
006b5f7f 425Int_t AliITStrackV2::Invariant() const {
426 //------------------------------------------------------------------
427 // This function is for debugging purpose only
428 //------------------------------------------------------------------
7f6ddf58 429 Int_t n=GetNumberOfClusters();
430
a9a2d814 431 if (TMath::Abs(fP2)>=0.9999){
8676d691 432 if (n>kWARN) Warning("Invariant","fP2=%f\n",fP2);
a9a2d814 433 return 0;
434 }
435 if (fC00<=0 || fC00>9.) {
8676d691 436 if (n>kWARN) Warning("Invariant","fC00=%f\n",fC00);
a9a2d814 437 return 0;
438 }
439 if (fC11<=0 || fC11>9.) {
8676d691 440 if (n>kWARN) Warning("Invariant","fC11=%f\n",fC11);
a9a2d814 441 return 0;
442 }
443 if (fC22<=0 || fC22>1.) {
8676d691 444 if (n>kWARN) Warning("Invariant","fC22=%f\n",fC22);
a9a2d814 445 return 0;
446 }
447 if (fC33<=0 || fC33>1.) {
8676d691 448 if (n>kWARN) Warning("Invariant","fC33=%f\n",fC33);
a9a2d814 449 return 0;
450 }
451 if (fC44<=0 || fC44>6e-5) {
8676d691 452 if (n>kWARN) Warning("Invariant","fC44=%f\n",fC44);
a9a2d814 453 return 0;
14825d5a 454 }
14825d5a 455 return 1;
006b5f7f 456}
457
458//____________________________________________________________________________
a9a2d814 459Int_t AliITStrackV2::Propagate(Double_t alp,Double_t xk) {
006b5f7f 460 //------------------------------------------------------------------
461 //This function propagates a track
462 //------------------------------------------------------------------
33329428 463 Double_t alpha=fAlpha, x=fX;
006b5f7f 464 Double_t p0=fP0,p1=fP1,p2=fP2,p3=fP3,p4=fP4;
465 Double_t c00=fC00;
466 Double_t c10=fC10, c11=fC11;
467 Double_t c20=fC20, c21=fC21, c22=fC22;
468 Double_t c30=fC30, c31=fC31, c32=fC32, c33=fC33;
469 Double_t c40=fC40, c41=fC41, c42=fC42, c43=fC43, c44=fC44;
470
33329428 471 if (alp < -TMath::Pi()) alp += 2*TMath::Pi();
472 else if (alp >= TMath::Pi()) alp -= 2*TMath::Pi();
473 Double_t ca=TMath::Cos(alp-fAlpha), sa=TMath::Sin(alp-fAlpha);
474 Double_t sf=fP2, cf=TMath::Sqrt(1.- fP2*fP2);
006b5f7f 475
33329428 476 // **** rotation **********************
477 {
006b5f7f 478 fAlpha = alp;
33329428 479 fX = x*ca + p0*sa;
480 fP0= -x*sa + p0*ca;
481 fP2= sf*ca - cf*sa;
482
c219fd63 483 Double_t rr=(ca+sf/cf*sa);
484
485 fC00 *= (ca*ca);
486 fC10 *= ca;
487 fC20 *= ca*rr;
488 fC30 *= ca;
489 fC40 *= ca;
490 //fC11 = fC11;
491 fC21 *= rr;
492 //fC31 = fC31;
493 //fC41 = fC41;
494 fC22 *= rr*rr;
495 fC32 *= rr;
496 fC42 *= rr;
497 //fC33=fC33;
498 //fC43=fC43;
499 //fC44=fC44;
500
33329428 501 }
006b5f7f 502
33329428 503 // **** translation ******************
504 {
505 Double_t dx=xk-fX;
006b5f7f 506 Double_t f1=fP2, f2=f1 + fP4*dx;
a9a2d814 507 if (TMath::Abs(f2) >= 0.9999) {
006b5f7f 508 Int_t n=GetNumberOfClusters();
7f6ddf58 509 if (n>kWARN)
8676d691 510 Warning("Propagate","Propagation failed (%d) !\n",n);
006b5f7f 511 return 0;
512 }
33329428 513 Double_t r1=TMath::Sqrt(1.- f1*f1), r2=TMath::Sqrt(1.- f2*f2);
006b5f7f 514
33329428 515 fX=xk;
006b5f7f 516 fP0 += dx*(f1+f2)/(r1+r2);
517 fP1 += dx*(f1+f2)/(f1*r2 + f2*r1)*fP3;
518 fP2 += dx*fP4;
519
c219fd63 520 //f = F - 1
521
522 Double_t f02= dx/(r1*r1*r1);
523 Double_t f04=0.5*dx*dx/(r1*r1*r1);
524 Double_t f12= dx*fP3*f1/(r1*r1*r1);
525 Double_t f14=0.5*dx*dx*fP3*f1/(r1*r1*r1);
526 Double_t f13= dx/r1;
527 Double_t f24= dx;
528
529 //b = C*ft
530 Double_t b00=f02*fC20 + f04*fC40, b01=f12*fC20 + f14*fC40 + f13*fC30;
531 Double_t b02=f24*fC40;
532 Double_t b10=f02*fC21 + f04*fC41, b11=f12*fC21 + f14*fC41 + f13*fC31;
533 Double_t b12=f24*fC41;
534 Double_t b20=f02*fC22 + f04*fC42, b21=f12*fC22 + f14*fC42 + f13*fC32;
535 Double_t b22=f24*fC42;
536 Double_t b40=f02*fC42 + f04*fC44, b41=f12*fC42 + f14*fC44 + f13*fC43;
537 Double_t b42=f24*fC44;
538 Double_t b30=f02*fC32 + f04*fC43, b31=f12*fC32 + f14*fC43 + f13*fC33;
539 Double_t b32=f24*fC43;
540
541 //a = f*b = f*C*ft
542 Double_t a00=f02*b20+f04*b40,a01=f02*b21+f04*b41,a02=f02*b22+f04*b42;
543 Double_t a11=f12*b21+f14*b41+f13*b31,a12=f12*b22+f14*b42+f13*b32;
544 Double_t a22=f24*b42;
545
546 //F*C*Ft = C + (b + bt + a)
547 fC00 += b00 + b00 + a00;
548 fC10 += b10 + b01 + a01;
549 fC20 += b20 + b02 + a02;
550 fC30 += b30;
551 fC40 += b40;
552 fC11 += b11 + b11 + a11;
553 fC21 += b21 + b12 + a12;
554 fC31 += b31;
555 fC41 += b41;
556 fC22 += b22 + b22 + a22;
557 fC32 += b32;
558 fC42 += b42;
006b5f7f 559
006b5f7f 560 if (!Invariant()) {
33329428 561 fAlpha=alpha;
562 fX=x;
006b5f7f 563 fP0=p0; fP1=p1; fP2=p2; fP3=p3; fP4=p4;
564 fC00=c00;
565 fC10=c10; fC11=c11;
566 fC20=c20; fC21=c21; fC22=c22;
567 fC30=c30; fC31=c31; fC32=c32; fC33=c33;
568 fC40=c40; fC41=c41; fC42=c42; fC43=c43; fC44=c44;
569 return 0;
570 }
33329428 571 }
006b5f7f 572
573 return 1;
574}
575
c219fd63 576
a6a6f799 577Double_t AliITStrackV2::GetD(Double_t x, Double_t y) const {
006b5f7f 578 //------------------------------------------------------------------
a6a6f799 579 // This function calculates the transverse impact parameter
580 // with respect to a point with global coordinates (x,y)
006b5f7f 581 //------------------------------------------------------------------
a6a6f799 582 Double_t xt=fX, yt=fP0;
583
584 Double_t sn=TMath::Sin(fAlpha), cs=TMath::Cos(fAlpha);
585 Double_t a = x*cs + y*sn;
586 y = -x*sn + y*cs; x=a;
587 xt-=x; yt-=y;
588
589 sn=fP4*xt - fP2; cs=fP4*yt + TMath::Sqrt(1.- fP2*fP2);
590 a=2*(xt*fP2 - yt*TMath::Sqrt(1.- fP2*fP2))-fP4*(xt*xt + yt*yt);
006b5f7f 591 if (fP4<0) a=-a;
592 return a/(1 + TMath::Sqrt(sn*sn + cs*cs));
a9a2d814 593}
006b5f7f 594
babd135a 595Double_t AliITStrackV2::GetZat(Double_t x) const {
596 //------------------------------------------------------------------
597 // This function calculates the z at given x point - in current coordinate system
598 //------------------------------------------------------------------
599 Double_t x1=fX, x2=x, dx=x2-x1;
600 //
601 Double_t f1=fP2, f2=f1 + fP4*dx;
602 if (TMath::Abs(f2) >= 0.9999) {
603 return 10000000;
604 }
babd135a 605 Double_t r1=sqrt(1.- f1*f1), r2=sqrt(1.- f2*f2);
606 Double_t z = fP1 + dx*(f1+f2)/(f1*r2 + f2*r1)*fP3;
607 return z;
608}
609
610
611
612
8676d691 613Int_t AliITStrackV2::Improve(Double_t x0,Double_t xyz[3],Double_t ers[3]) {
a9a2d814 614 //------------------------------------------------------------------
615 //This function improves angular track parameters
616 //------------------------------------------------------------------
8676d691 617 Double_t cs=TMath::Cos(fAlpha), sn=TMath::Sin(fAlpha);
618 //Double_t xv = xyz[0]*cs + xyz[1]*sn; // vertex
619 Double_t yv =-xyz[0]*sn + xyz[1]*cs; // in the
620 Double_t zv = xyz[2]; // local frame
a9a2d814 621 Double_t dy=fP0-yv, dz=fP1-zv;
622 Double_t r2=fX*fX+dy*dy;
623 Double_t p2=(1.+ GetTgl()*GetTgl())/(Get1Pt()*Get1Pt());
624 Double_t beta2=p2/(p2 + GetMass()*GetMass());
625 x0*=TMath::Sqrt((1.+ GetTgl()*GetTgl())/(1.- GetSnp()*GetSnp()));
8676d691 626 Double_t theta2=14.1*14.1/(beta2*p2*1e6)*x0;
627 //Double_t theta2=1.0259e-6*14*14/28/(beta2*p2)*x0*9.36*2.33;
a9a2d814 628 {
25d954aa 629 Double_t dummy=4/r2-fP4*fP4;
630 if (dummy < 0) return 0;
631 Double_t parp=0.5*(fP4*fX + dy*TMath::Sqrt(dummy));
a9a2d814 632 Double_t sigma2p = theta2*(1.- GetSnp()*GetSnp())*(1. + GetTgl()*GetTgl());
633 sigma2p += fC00/r2*(1.- dy*dy/r2)*(1.- dy*dy/r2);
8676d691 634 sigma2p += ers[1]*ers[1]/r2;
a9a2d814 635 sigma2p += 0.25*fC44*fX*fX;
636 Double_t eps2p=sigma2p/(fC22+sigma2p);
637 fP0 += fC20/(fC22+sigma2p)*(parp-fP2);
638 fP2 = eps2p*fP2 + (1-eps2p)*parp;
639 fC22 *= eps2p;
640 fC20 *= eps2p;
641 }
642 {
643 Double_t parl=0.5*fP4*dz/TMath::ASin(0.5*fP4*TMath::Sqrt(r2));
644 Double_t sigma2l=theta2;
645 sigma2l += fC11/r2+fC00*dy*dy*dz*dz/(r2*r2*r2);
8676d691 646 sigma2l += ers[2]*ers[2]/r2;
a9a2d814 647 Double_t eps2l=sigma2l/(fC33+sigma2l);
648 fP1 += fC31/(fC33+sigma2l)*(parl-fP3);
649 fP4 += fC43/(fC33+sigma2l)*(parl-fP3);
650 fP3 = eps2l*fP3 + (1-eps2l)*parl;
651 fC33 *= eps2l; fC43 *= eps2l;
652 fC31 *= eps2l;
653 }
654 if (!Invariant()) return 0;
655 return 1;
656}
006b5f7f 657
14825d5a 658void AliITStrackV2::ResetCovariance() {
659 //------------------------------------------------------------------
660 //This function makes a track forget its history :)
661 //------------------------------------------------------------------
662
663 fC00*=10.;
664 fC10=0.; fC11*=10.;
665 fC20=0.; fC21=0.; fC22*=10.;
666 fC30=0.; fC31=0.; fC32=0.; fC33*=10.;
667 fC40=0.; fC41=0.; fC42=0.; fC43=0.; fC44*=10.;
668
669}
23efe5f1 670
671void AliITStrackV2::CookdEdx(Double_t low, Double_t up) {
672 //-----------------------------------------------------------------
a9a2d814 673 // This function calculates dE/dX within the "low" and "up" cuts.
674 // Origin: Boris Batyunya, JINR, Boris.Batiounia@cern.ch
23efe5f1 675 //-----------------------------------------------------------------
23efe5f1 676 // The clusters order is: SSD-2, SSD-1, SDD-2, SDD-1, SPD-2, SPD-1
c0dd5278 677
678 Int_t i;
679 Int_t nc=0;
680 for (i=0; i<GetNumberOfClusters(); i++) {
681 Int_t idx=GetClusterIndex(i);
682 idx=(idx&0xf0000000)>>28;
683 if (idx>1) nc++; // Take only SSD and SDD
684 }
23efe5f1 685
a9a2d814 686 Int_t swap;//stupid sorting
23efe5f1 687 do {
688 swap=0;
689 for (i=0; i<nc-1; i++) {
690 if (fdEdxSample[i]<=fdEdxSample[i+1]) continue;
691 Float_t tmp=fdEdxSample[i];
692 fdEdxSample[i]=fdEdxSample[i+1]; fdEdxSample[i+1]=tmp;
693 swap++;
694 }
695 } while (swap);
696
a9a2d814 697 Int_t nl=Int_t(low*nc), nu=Int_t(up*nc); //b.b. to take two lowest dEdX
698 // values from four ones choose
699 // nu=2
23efe5f1 700 Float_t dedx=0;
a9a2d814 701 for (i=nl; i<nu; i++) dedx += fdEdxSample[i];
c0dd5278 702 if (nu-nl>0) dedx /= (nu-nl);
23efe5f1 703
23efe5f1 704 SetdEdx(dedx);
705}