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