]> 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
PMD new raw data writer
[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
66e811e2 16///////////////////////////////////////////////////////////////////////////
006b5f7f 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
66e811e2 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{
ef7253ac 59 for(Int_t i=0; i<2*kMaxLayer; i++) fIndex[i]=-1;
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];
c0b978f0 84 if (c) t.GetConstrainedExternalParameters(fAlpha,x,p);
67c3dcbe 85 else t.GetExternalParameters(x,p);
c84a5e9e 86 fX=x;
006b5f7f 87 fP0=p[0];
c84a5e9e 88 fP1=p[1]; SaveLocalConvConst();
006b5f7f 89 fP2=p[2];
c84a5e9e 90 fP3=p[3]; x=GetLocalConvConst();
006b5f7f 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;
006b5f7f 112}
113
15dd636f 114void AliITStrackV2::UpdateESDtrack(ULong_t flags) const {
83d73500 115 fESDtrack->UpdateTrackParams(this,flags);
116}
117
006b5f7f 118//____________________________________________________________________________
119AliITStrackV2::AliITStrackV2(const AliITStrackV2& t) : AliKalmanTrack(t) {
120 //------------------------------------------------------------------
121 //Copy constructor
122 //------------------------------------------------------------------
123 fX=t.fX;
124 fAlpha=t.fAlpha;
125 fdEdx=t.fdEdx;
126
127 fP0=t.fP0; fP1=t.fP1; fP2=t.fP2; fP3=t.fP3; fP4=t.fP4;
128
129 fC00=t.fC00;
130 fC10=t.fC10; fC11=t.fC11;
131 fC20=t.fC20; fC21=t.fC21; fC22=t.fC22;
132 fC30=t.fC30; fC31=t.fC31; fC32=t.fC32; fC33=t.fC33;
133 fC40=t.fC40; fC41=t.fC41; fC42=t.fC42; fC43=t.fC43; fC44=t.fC44;
134
ef7253ac 135 Int_t i;
136 for (i=0; i<2*kMaxLayer; i++) fIndex[i]=t.fIndex[i];
137 for (i=0; i<4; i++) fdEdxSample[i]=t.fdEdxSample[i];
138
83d73500 139 fESDtrack=t.fESDtrack;
006b5f7f 140}
a9a2d814 141
006b5f7f 142//_____________________________________________________________________________
143Int_t AliITStrackV2::Compare(const TObject *o) const {
144 //-----------------------------------------------------------------
145 // This function compares tracks according to the their curvature
146 //-----------------------------------------------------------------
7f6ddf58 147 AliITStrackV2 *t=(AliITStrackV2*)o;
a9a2d814 148 //Double_t co=TMath::Abs(t->Get1Pt());
149 //Double_t c =TMath::Abs(Get1Pt());
15dd636f 150 Double_t co=t->GetSigmaY2()*t->GetSigmaZ2();
151 Double_t c =GetSigmaY2()*GetSigmaZ2();
006b5f7f 152 if (c>co) return 1;
153 else if (c<co) return -1;
154 return 0;
155}
156
157//_____________________________________________________________________________
158void AliITStrackV2::GetExternalCovariance(Double_t cc[15]) const {
159 //-------------------------------------------------------------------------
160 // This function returns an external representation of the covriance matrix.
161 // (See comments in AliTPCtrack.h about external track representation)
162 //-------------------------------------------------------------------------
c84a5e9e 163 Double_t a=GetLocalConvConst();
006b5f7f 164
165 cc[0 ]=fC00;
166 cc[1 ]=fC10; cc[2 ]=fC11;
167 cc[3 ]=fC20; cc[4 ]=fC21; cc[5 ]=fC22;
168 cc[6 ]=fC30; cc[7 ]=fC31; cc[8 ]=fC32; cc[9 ]=fC33;
169 cc[10]=fC40*a; cc[11]=fC41*a; cc[12]=fC42*a; cc[13]=fC43*a; cc[14]=fC44*a*a;
170}
171
172//____________________________________________________________________________
a9a2d814 173Int_t AliITStrackV2::PropagateToVertex(Double_t d,Double_t x0) {
006b5f7f 174 //------------------------------------------------------------------
175 //This function propagates a track to the minimal distance from the origin
176 //------------------------------------------------------------------
5b316485 177 //Double_t xv=fP2*(fX*fP2 - fP0*TMath::Sqrt(1.- fP2*fP2)); //linear approxim.
178 Double_t tgf=-(fP4*fX - fP2)/(fP4*fP0 + TMath::Sqrt(1 - fP2*fP2));
179 Double_t snf=tgf/TMath::Sqrt(1.+ tgf*tgf);
180 Double_t xv=(snf - fP2)/fP4 + fX;
181 return PropagateTo(xv,d,x0);
006b5f7f 182}
183
184//____________________________________________________________________________
185Int_t AliITStrackV2::
186GetGlobalXYZat(Double_t xk, Double_t &x, Double_t &y, Double_t &z) const {
187 //------------------------------------------------------------------
188 //This function returns a track position in the global system
189 //------------------------------------------------------------------
190 Double_t dx=xk-fX;
191 Double_t f1=fP2, f2=f1 + fP4*dx;
a9a2d814 192 if (TMath::Abs(f2) >= 0.9999) {
7f6ddf58 193 Int_t n=GetNumberOfClusters();
194 if (n>kWARN)
8676d691 195 Warning("GetGlobalXYZat","Propagation failed (%d) !\n",n);
006b5f7f 196 return 0;
197 }
198
199 Double_t r1=sqrt(1.- f1*f1), r2=sqrt(1.- f2*f2);
200
201 Double_t yk = fP0 + dx*(f1+f2)/(r1+r2);
202 Double_t zk = fP1 + dx*(f1+f2)/(f1*r2 + f2*r1)*fP3;
203
204 Double_t cs=TMath::Cos(fAlpha), sn=TMath::Sin(fAlpha);
205 x = xk*cs - yk*sn;
206 y = xk*sn + yk*cs;
207 z = zk;
208
209 return 1;
210}
211
212//_____________________________________________________________________________
213Double_t AliITStrackV2::GetPredictedChi2(const AliCluster *c) const
214{
215 //-----------------------------------------------------------------
216 // This function calculates a predicted chi2 increment.
217 //-----------------------------------------------------------------
218 Double_t r00=c->GetSigmaY2(), r01=0., r11=c->GetSigmaZ2();
219 r00+=fC00; r01+=fC10; r11+=fC11;
babd135a 220 //
006b5f7f 221 Double_t det=r00*r11 - r01*r01;
7f6ddf58 222 if (TMath::Abs(det) < 1.e-30) {
006b5f7f 223 Int_t n=GetNumberOfClusters();
7f6ddf58 224 if (n>kWARN)
8676d691 225 Warning("GetPredictedChi2","Singular matrix (%d) !\n",n);
006b5f7f 226 return 1e10;
227 }
228 Double_t tmp=r00; r00=r11; r11=tmp; r01=-r01;
33329428 229
006b5f7f 230 Double_t dy=c->GetY() - fP0, dz=c->GetZ() - fP1;
33329428 231
006b5f7f 232 return (dy*r00*dy + 2*r01*dy*dz + dz*r11*dz)/det;
233}
234
006b5f7f 235//____________________________________________________________________________
a9a2d814 236Int_t AliITStrackV2::CorrectForMaterial(Double_t d, Double_t x0) {
237 //------------------------------------------------------------------
238 //This function corrects the track parameters for crossed material
239 //------------------------------------------------------------------
e43c066c 240 Double_t p2=(1.+ fP3*fP3)/(Get1Pt()*Get1Pt());
15dd636f 241 Double_t beta2=p2/(p2 + GetMass()*GetMass());
a9a2d814 242 d*=TMath::Sqrt((1.+ fP3*fP3)/(1.- fP2*fP2));
243
244 //Multiple scattering******************
245 if (d!=0) {
33329428 246 Double_t theta2=14.1*14.1/(beta2*p2*1e6)*TMath::Abs(d);
247 //Double_t theta2=1.0259e-6*14*14/28/(beta2*p2)*TMath::Abs(d)*9.36*2.33;
a9a2d814 248 fC22 += theta2*(1.- fP2*fP2)*(1. + fP3*fP3);
249 fC33 += theta2*(1. + fP3*fP3)*(1. + fP3*fP3);
250 fC43 += theta2*fP3*fP4*(1. + fP3*fP3);
251 fC44 += theta2*fP3*fP4*fP3*fP4;
252 }
253
254 //Energy losses************************
255 if (x0!=0.) {
256 d*=x0;
15dd636f 257 Double_t dE=0.153e-3/beta2*(log(5940*beta2/(1-beta2)) - beta2)*d;
258 if (beta2/(1-beta2)>3.5*3.5)
83d73500 259 dE=0.153e-3/beta2*(log(3.5*5940)+0.5*log(beta2/(1-beta2)) - beta2)*d;
15dd636f 260
261 fP4*=(1.- TMath::Sqrt(p2+GetMass()*GetMass())/p2*dE);
a9a2d814 262 }
263
264 if (!Invariant()) return 0;
265
266 return 1;
267}
268
269//____________________________________________________________________________
270Int_t AliITStrackV2::PropagateTo(Double_t xk, Double_t d, Double_t x0) {
006b5f7f 271 //------------------------------------------------------------------
272 //This function propagates a track
273 //------------------------------------------------------------------
274 Double_t x1=fX, x2=xk, dx=x2-x1;
275 Double_t f1=fP2, f2=f1 + fP4*dx;
c84a5e9e 276 if (TMath::Abs(f2) >= 0.98) {
81e97e0d 277 // MI change - don't propagate highly inclined tracks
278 // covariance matrix distorted
c84a5e9e 279 //Int_t n=GetNumberOfClusters();
280 //if (n>kWARN)
281 // Warning("PropagateTo","Propagation failed !\n",n);
006b5f7f 282 return 0;
283 }
c84a5e9e 284 Double_t lcc=GetLocalConvConst();
006b5f7f 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
c84a5e9e 337 //Change of the magnetic field *************
338 SaveLocalConvConst();
339 fP4*=lcc/GetLocalConvConst();
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;
c84a5e9e 509 if (TMath::Abs(f2) >= 0.98) {
81e97e0d 510 // don't propagate highly inclined tracks MI
006b5f7f 511 return 0;
512 }
c84a5e9e 513 // Int_t n=GetNumberOfClusters();
514 // if (n>kWARN)
515 // Warning("Propagate","Propagation failed (%d) !\n",n);
516 // return 0;
517 //}
518 Double_t lcc=GetLocalConvConst();
519
33329428 520 Double_t r1=TMath::Sqrt(1.- f1*f1), r2=TMath::Sqrt(1.- f2*f2);
c84a5e9e 521
33329428 522 fX=xk;
006b5f7f 523 fP0 += dx*(f1+f2)/(r1+r2);
524 fP1 += dx*(f1+f2)/(f1*r2 + f2*r1)*fP3;
525 fP2 += dx*fP4;
526
c84a5e9e 527 //Change of the magnetic field *************
528 SaveLocalConvConst();
529 fP4*=lcc/GetLocalConvConst();
530
c219fd63 531 //f = F - 1
532
533 Double_t f02= dx/(r1*r1*r1);
534 Double_t f04=0.5*dx*dx/(r1*r1*r1);
535 Double_t f12= dx*fP3*f1/(r1*r1*r1);
536 Double_t f14=0.5*dx*dx*fP3*f1/(r1*r1*r1);
537 Double_t f13= dx/r1;
538 Double_t f24= dx;
539
540 //b = C*ft
541 Double_t b00=f02*fC20 + f04*fC40, b01=f12*fC20 + f14*fC40 + f13*fC30;
542 Double_t b02=f24*fC40;
543 Double_t b10=f02*fC21 + f04*fC41, b11=f12*fC21 + f14*fC41 + f13*fC31;
544 Double_t b12=f24*fC41;
545 Double_t b20=f02*fC22 + f04*fC42, b21=f12*fC22 + f14*fC42 + f13*fC32;
546 Double_t b22=f24*fC42;
547 Double_t b40=f02*fC42 + f04*fC44, b41=f12*fC42 + f14*fC44 + f13*fC43;
548 Double_t b42=f24*fC44;
549 Double_t b30=f02*fC32 + f04*fC43, b31=f12*fC32 + f14*fC43 + f13*fC33;
550 Double_t b32=f24*fC43;
551
552 //a = f*b = f*C*ft
553 Double_t a00=f02*b20+f04*b40,a01=f02*b21+f04*b41,a02=f02*b22+f04*b42;
554 Double_t a11=f12*b21+f14*b41+f13*b31,a12=f12*b22+f14*b42+f13*b32;
555 Double_t a22=f24*b42;
556
557 //F*C*Ft = C + (b + bt + a)
558 fC00 += b00 + b00 + a00;
559 fC10 += b10 + b01 + a01;
560 fC20 += b20 + b02 + a02;
561 fC30 += b30;
562 fC40 += b40;
563 fC11 += b11 + b11 + a11;
564 fC21 += b21 + b12 + a12;
565 fC31 += b31;
566 fC41 += b41;
567 fC22 += b22 + b22 + a22;
568 fC32 += b32;
569 fC42 += b42;
006b5f7f 570
006b5f7f 571 if (!Invariant()) {
33329428 572 fAlpha=alpha;
573 fX=x;
006b5f7f 574 fP0=p0; fP1=p1; fP2=p2; fP3=p3; fP4=p4;
575 fC00=c00;
576 fC10=c10; fC11=c11;
577 fC20=c20; fC21=c21; fC22=c22;
578 fC30=c30; fC31=c31; fC32=c32; fC33=c33;
579 fC40=c40; fC41=c41; fC42=c42; fC43=c43; fC44=c44;
580 return 0;
581 }
33329428 582 }
006b5f7f 583
584 return 1;
585}
586
c219fd63 587
a6a6f799 588Double_t AliITStrackV2::GetD(Double_t x, Double_t y) const {
006b5f7f 589 //------------------------------------------------------------------
a6a6f799 590 // This function calculates the transverse impact parameter
591 // with respect to a point with global coordinates (x,y)
006b5f7f 592 //------------------------------------------------------------------
a6a6f799 593 Double_t xt=fX, yt=fP0;
594
595 Double_t sn=TMath::Sin(fAlpha), cs=TMath::Cos(fAlpha);
596 Double_t a = x*cs + y*sn;
597 y = -x*sn + y*cs; x=a;
598 xt-=x; yt-=y;
599
600 sn=fP4*xt - fP2; cs=fP4*yt + TMath::Sqrt(1.- fP2*fP2);
601 a=2*(xt*fP2 - yt*TMath::Sqrt(1.- fP2*fP2))-fP4*(xt*xt + yt*yt);
791f9a2a 602 return -a/(1 + TMath::Sqrt(sn*sn + cs*cs));
603}
604
605void
606AliITStrackV2::GetDZ(Double_t x, Double_t y, Double_t z, Float_t dz[2]) const {
607 //------------------------------------------------------------------
608 // This function calculates the transverse and longitudinal impact parameters
609 // with respect to a point with global coordinates (x,y,z)
610 //------------------------------------------------------------------
611 Double_t f1 = fP2, r1 = TMath::Sqrt(1. - f1*f1);
612 Double_t xt=fX, yt=fP0;
613 Double_t sn=TMath::Sin(fAlpha), cs=TMath::Cos(fAlpha);
614 Double_t a = x*cs + y*sn;
615 y = -x*sn + y*cs; x=a;
616 xt-=x; yt-=y;
617
618 Double_t rp4=fP4;
619 if (TMath::Abs(rp4) < kAlmost0) {
620 dz[0] = -(xt*f1 - yt*r1);
621 dz[1] = fP1 + (dz[0]*f1 - xt)/r1*fP3 - z;
622 return;
623 }
624
625 sn=rp4*xt - f1; cs=rp4*yt + r1;
626 a=2*(xt*f1 - yt*r1)-rp4*(xt*xt + yt*yt);
627 Double_t rr=TMath::Sqrt(sn*sn + cs*cs);
628 dz[0] = -a/(1 + rr);
629 Double_t f2 = -sn/rr, r2 = TMath::Sqrt(1. - f2*f2);
630 dz[1] = fP1 + fP3/rp4*TMath::ASin(f2*r1 - f1*r2) - z;
a9a2d814 631}
006b5f7f 632
babd135a 633Double_t AliITStrackV2::GetZat(Double_t x) const {
634 //------------------------------------------------------------------
635 // This function calculates the z at given x point - in current coordinate system
636 //------------------------------------------------------------------
637 Double_t x1=fX, x2=x, dx=x2-x1;
638 //
639 Double_t f1=fP2, f2=f1 + fP4*dx;
640 if (TMath::Abs(f2) >= 0.9999) {
641 return 10000000;
642 }
babd135a 643 Double_t r1=sqrt(1.- f1*f1), r2=sqrt(1.- f2*f2);
644 Double_t z = fP1 + dx*(f1+f2)/(f1*r2 + f2*r1)*fP3;
645 return z;
646}
647
648
649
650
8676d691 651Int_t AliITStrackV2::Improve(Double_t x0,Double_t xyz[3],Double_t ers[3]) {
a9a2d814 652 //------------------------------------------------------------------
653 //This function improves angular track parameters
654 //------------------------------------------------------------------
8676d691 655 Double_t cs=TMath::Cos(fAlpha), sn=TMath::Sin(fAlpha);
656 //Double_t xv = xyz[0]*cs + xyz[1]*sn; // vertex
657 Double_t yv =-xyz[0]*sn + xyz[1]*cs; // in the
658 Double_t zv = xyz[2]; // local frame
a9a2d814 659 Double_t dy=fP0-yv, dz=fP1-zv;
660 Double_t r2=fX*fX+dy*dy;
661 Double_t p2=(1.+ GetTgl()*GetTgl())/(Get1Pt()*Get1Pt());
662 Double_t beta2=p2/(p2 + GetMass()*GetMass());
663 x0*=TMath::Sqrt((1.+ GetTgl()*GetTgl())/(1.- GetSnp()*GetSnp()));
8676d691 664 Double_t theta2=14.1*14.1/(beta2*p2*1e6)*x0;
665 //Double_t theta2=1.0259e-6*14*14/28/(beta2*p2)*x0*9.36*2.33;
a9a2d814 666 {
25d954aa 667 Double_t dummy=4/r2-fP4*fP4;
668 if (dummy < 0) return 0;
669 Double_t parp=0.5*(fP4*fX + dy*TMath::Sqrt(dummy));
a9a2d814 670 Double_t sigma2p = theta2*(1.- GetSnp()*GetSnp())*(1. + GetTgl()*GetTgl());
671 sigma2p += fC00/r2*(1.- dy*dy/r2)*(1.- dy*dy/r2);
8676d691 672 sigma2p += ers[1]*ers[1]/r2;
a9a2d814 673 sigma2p += 0.25*fC44*fX*fX;
674 Double_t eps2p=sigma2p/(fC22+sigma2p);
675 fP0 += fC20/(fC22+sigma2p)*(parp-fP2);
676 fP2 = eps2p*fP2 + (1-eps2p)*parp;
677 fC22 *= eps2p;
678 fC20 *= eps2p;
679 }
680 {
681 Double_t parl=0.5*fP4*dz/TMath::ASin(0.5*fP4*TMath::Sqrt(r2));
682 Double_t sigma2l=theta2;
683 sigma2l += fC11/r2+fC00*dy*dy*dz*dz/(r2*r2*r2);
8676d691 684 sigma2l += ers[2]*ers[2]/r2;
a9a2d814 685 Double_t eps2l=sigma2l/(fC33+sigma2l);
686 fP1 += fC31/(fC33+sigma2l)*(parl-fP3);
687 fP4 += fC43/(fC33+sigma2l)*(parl-fP3);
688 fP3 = eps2l*fP3 + (1-eps2l)*parl;
689 fC33 *= eps2l; fC43 *= eps2l;
690 fC31 *= eps2l;
691 }
692 if (!Invariant()) return 0;
693 return 1;
694}
006b5f7f 695
14825d5a 696void AliITStrackV2::ResetCovariance() {
697 //------------------------------------------------------------------
698 //This function makes a track forget its history :)
699 //------------------------------------------------------------------
700
701 fC00*=10.;
702 fC10=0.; fC11*=10.;
703 fC20=0.; fC21=0.; fC22*=10.;
704 fC30=0.; fC31=0.; fC32=0.; fC33*=10.;
705 fC40=0.; fC41=0.; fC42=0.; fC43=0.; fC44*=10.;
706
707}
23efe5f1 708
709void AliITStrackV2::CookdEdx(Double_t low, Double_t up) {
710 //-----------------------------------------------------------------
a9a2d814 711 // This function calculates dE/dX within the "low" and "up" cuts.
712 // Origin: Boris Batyunya, JINR, Boris.Batiounia@cern.ch
23efe5f1 713 //-----------------------------------------------------------------
23efe5f1 714 // The clusters order is: SSD-2, SSD-1, SDD-2, SDD-1, SPD-2, SPD-1
c0dd5278 715
716 Int_t i;
717 Int_t nc=0;
718 for (i=0; i<GetNumberOfClusters(); i++) {
719 Int_t idx=GetClusterIndex(i);
720 idx=(idx&0xf0000000)>>28;
721 if (idx>1) nc++; // Take only SSD and SDD
722 }
23efe5f1 723
a9a2d814 724 Int_t swap;//stupid sorting
23efe5f1 725 do {
726 swap=0;
727 for (i=0; i<nc-1; i++) {
728 if (fdEdxSample[i]<=fdEdxSample[i+1]) continue;
729 Float_t tmp=fdEdxSample[i];
730 fdEdxSample[i]=fdEdxSample[i+1]; fdEdxSample[i+1]=tmp;
731 swap++;
732 }
733 } while (swap);
734
a9a2d814 735 Int_t nl=Int_t(low*nc), nu=Int_t(up*nc); //b.b. to take two lowest dEdX
736 // values from four ones choose
737 // nu=2
23efe5f1 738 Float_t dedx=0;
a9a2d814 739 for (i=nl; i<nu; i++) dedx += fdEdxSample[i];
c0dd5278 740 if (nu-nl>0) dedx /= (nu-nl);
23efe5f1 741
23efe5f1 742 SetdEdx(dedx);
743}
c7bafca9 744
745Double_t AliITStrackV2::
746PropagateToDCA(AliKalmanTrack *p, Double_t d, Double_t x0) {
747 //--------------------------------------------------------------
748 // Propagates this track and the argument track to the position of the
749 // distance of closest approach.
750 // Returns the (weighed !) distance of closest approach.
751 //--------------------------------------------------------------
752 Double_t xthis, xp, dca;
753 {
754 //Temporary solution
755 Double_t b=1./GetLocalConvConst()/kB2C;
756 AliExternalTrackParam dummy1(*this), dummy2(*p);
757 dca=dummy1.GetDCA(&dummy2,b,xthis,xp);
758 }
759 if (!PropagateTo(xthis,d,x0)) {
760 //AliWarning(" propagation failed !");
761 return 1e+33;
762 }
763
764 if (!p->PropagateTo(xp,d,x0)) {
765 //AliWarning(" propagation failed !";
766 return 1e+33;
767 }
768
769 return dca;
770}
5773defd 771
772Double_t AliITStrackV2::Get1Pt() const {
773 //--------------------------------------------------------------
774 // Returns the inverse Pt (1/GeV/c)
775 // (or 1/"most probable pt", if the field is too weak)
776 //--------------------------------------------------------------
777 if (TMath::Abs(GetLocalConvConst()) > kVeryBigConvConst)
778 return 1./kMostProbableMomentum/TMath::Sqrt(1.+ GetTgl()*GetTgl());
779 return (TMath::Sign(1e-9,fP4) + fP4)*GetLocalConvConst();
780}