]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TRD/AliTRDtrack.cxx
Effective C++ warnings
[u/mrichter/AliRoot.git] / TRD / AliTRDtrack.cxx
CommitLineData
46d29e70 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
0fa7dfa7 16/* $Id$ */
46d29e70 17
a2cb5b3d 18#include <Riostream.h>
53c17fbf 19#include <TMath.h>
e1e6896f 20#include <TVector2.h>
46d29e70 21
53c17fbf 22#include "AliESDtrack.h"
46d29e70 23#include "AliTRDgeometry.h"
24#include "AliTRDcluster.h"
25#include "AliTRDtrack.h"
53c17fbf 26#include "AliTRDtracklet.h"
b3a5a838 27
46d29e70 28ClassImp(AliTRDtrack)
29
53c17fbf 30///////////////////////////////////////////////////////////////////////////////
31// //
32// Represents a reconstructed TRD track //
33// Local TRD Kalman track //
34// //
35///////////////////////////////////////////////////////////////////////////////
7ad19338 36
15ed8ba1 37//_____________________________________________________________________________
38AliTRDtrack::AliTRDtrack()
39 :AliKalmanTrack()
40 ,fSeedLab(-1)
41 ,fdEdx(0)
15ed8ba1 42 ,fDE(0)
43 ,fAlpha(0)
44 ,fX(0)
45 ,fStopped(kFALSE)
46 ,fY(0)
47 ,fZ(0)
48 ,fE(0)
49 ,fT(0)
50 ,fC(0)
51 ,fCyy(1e10)
52 ,fCzy(0)
53 ,fCzz(1e10)
54 ,fCey(0)
55 ,fCez(0)
56 ,fCee(1e10)
57 ,fCty(0)
58 ,fCtz(0)
59 ,fCte(0)
60 ,fCtt(1e10)
61 ,fCcy(0)
62 ,fCcz(0)
63 ,fCce(0)
64 ,fCct(0)
65 ,fCcc(1e10)
66 ,fLhElectron(0)
67 ,fNWrong(0)
68 ,fNRotate(0)
69 ,fNCross(0)
70 ,fNExpected(0)
71 ,fNLast(0)
72 ,fNExpectedLast(0)
73 ,fNdedx(0)
74 ,fChi2Last(1e10)
75 ,fBackupTrack(0x0)
f146da82 76{
15ed8ba1 77 //
78 // AliTRDtrack default constructor
79 //
80
81 Int_t i = 0;
82 Int_t j = 0;
83 UInt_t k = 0;
84
85 for (i = 0; i < kNplane; i++) {
86 for (j = 0; j < kNslice; j++) {
6d45eaef 87 fdEdxPlane[i][j] = 0;
88 }
f146da82 89 fTimBinPlane[i] = -1;
90 }
15ed8ba1 91
92 for (k = 0; k < kMAXCLUSTERSPERTRACK; k++) {
93 fIndex[k] = 0;
94 fIndexBackup[k] = 0;
95 fdQdl[k] = 0;
96 }
97
98 for (i = 0; i < 3; i++) {
99 fBudget[i] = 0;
f146da82 100 }
15ed8ba1 101
f146da82 102}
6d45eaef 103
46d29e70 104//_____________________________________________________________________________
a819a5f7 105AliTRDtrack::AliTRDtrack(const AliTRDcluster *c, UInt_t index,
5443e65e 106 const Double_t xx[5], const Double_t cc[15],
15ed8ba1 107 Double_t xref, Double_t alpha)
108 :AliKalmanTrack()
109 ,fSeedLab(-1)
110 ,fdEdx(0.0)
15ed8ba1 111 ,fDE(0.0)
112 ,fAlpha(alpha)
113 ,fX(xref)
114 ,fStopped(kFALSE)
115 ,fY(xx[0])
116 ,fZ(xx[1])
117 ,fE(xx[2])
118 ,fT(xx[3])
119 ,fC(xx[4])
120 ,fCyy(cc[0])
121 ,fCzy(cc[1])
122 ,fCzz(cc[2])
123 ,fCey(cc[3])
124 ,fCez(cc[4])
125 ,fCee(cc[5])
126 ,fCty(cc[6])
127 ,fCtz(cc[7])
128 ,fCte(cc[8])
129 ,fCtt(cc[9])
130 ,fCcy(cc[10])
131 ,fCcz(cc[11])
132 ,fCce(cc[12])
133 ,fCct(cc[13])
134 ,fCcc(cc[14])
135 ,fLhElectron(0.0)
136 ,fNWrong(0)
137 ,fNRotate(0)
138 ,fNCross(0)
139 ,fNExpected(0)
140 ,fNLast(0)
141 ,fNExpectedLast(0)
142 ,fNdedx(0)
143 ,fChi2Last(1e10)
144 ,fBackupTrack(0x0)
145{
146 //
147 // AliTRDtrack main constructor
148 //
5443e65e 149
15ed8ba1 150 Int_t i = 0;
151 Int_t j = 0;
152 UInt_t k = 0;
46d29e70 153
d49e463b 154 if (fAlpha < -TMath::Pi()) {
155 fAlpha += 2.0 * TMath::Pi();
156 }
157 if (fAlpha >= TMath::Pi()) {
158 fAlpha -= 2.0 * TMath::Pi();
159 }
160
c84a5e9e 161 SaveLocalConvConst();
46d29e70 162
15ed8ba1 163 fIndex[0] = index;
5443e65e 164 SetNumberOfClusters(1);
165
15ed8ba1 166 for (i = 0; i < kNplane; i++) {
167 for (j = 0; j < kNslice; j++) {
6d45eaef 168 fdEdxPlane[i][j] = 0;
169 }
170 fTimBinPlane[i] = -1;
eab5961e 171 }
a2b90f83 172
5443e65e 173 Double_t q = TMath::Abs(c->GetQ());
15ed8ba1 174 Double_t s = fX * fC - fE;
175 Double_t t = fT;
176 if (s*s < 1.0) {
177 q *= TMath::Sqrt((1-s*s)/(1+t*t));
178 }
5443e65e 179 fdQdl[0] = q;
0d5b5c27 180
15ed8ba1 181 for (k = 1; k < kMAXCLUSTERSPERTRACK; k++) {
182 fdQdl[k] = 0;
183 fIndex[k] = 0;
184 fIndexBackup[k] = 0;
0d5b5c27 185 }
53c17fbf 186
15ed8ba1 187 for (i = 0; i < 3; i++) {
188 fBudget[i] = 0;
189 }
53c17fbf 190
46d29e70 191}
192
193//_____________________________________________________________________________
15ed8ba1 194AliTRDtrack::AliTRDtrack(const AliTRDtrack &t)
195 :AliKalmanTrack(t)
196 ,fSeedLab(t.fSeedLab)
197 ,fdEdx(t.fdEdx)
15ed8ba1 198 ,fDE(t.fDE)
199 ,fAlpha(t.fAlpha)
200 ,fX(t.fX)
201 ,fStopped(t.fStopped)
202 ,fY(t.fY)
203 ,fZ(t.fZ)
204 ,fE(t.fE)
205 ,fT(t.fT)
206 ,fC(t.fC)
207 ,fCyy(t.fCyy)
208 ,fCzy(t.fCzy)
209 ,fCzz(t.fCzz)
210 ,fCey(t.fCey)
211 ,fCez(t.fCez)
212 ,fCee(t.fCee)
213 ,fCty(t.fCty)
214 ,fCtz(t.fCtz)
215 ,fCte(t.fCte)
216 ,fCtt(t.fCtt)
217 ,fCcy(t.fCcy)
218 ,fCcz(t.fCcz)
219 ,fCce(t.fCce)
220 ,fCct(t.fCct)
221 ,fCcc(t.fCcc)
222 ,fLhElectron(0.0)
223 ,fNWrong(t.fNWrong)
224 ,fNRotate(t.fNRotate)
225 ,fNCross(t.fNCross)
226 ,fNExpected(t.fNExpected)
227 ,fNLast(t.fNLast)
228 ,fNExpectedLast(t.fNExpectedLast)
229 ,fNdedx(t.fNdedx)
230 ,fChi2Last(t.fChi2Last)
231 ,fBackupTrack(0x0)
53c17fbf 232{
46d29e70 233 //
234 // Copy constructor.
235 //
15ed8ba1 236
237 Int_t i = 0;
238 Int_t j = 0;
239 UInt_t k = 0;
240
241 for (i = 0; i < kNplane; i++) {
242 for (j = 0; j < kNslice; j++) {
6d45eaef 243 fdEdxPlane[i][j] = t.fdEdxPlane[i][j];
244 }
245 fTimBinPlane[i] = t.fTimBinPlane[i];
246 fTracklets[i] = t.fTracklets[i];
eab5961e 247 }
46d29e70 248
15ed8ba1 249 Int_t n = t.GetNumberOfClusters();
250 for (i = 0; i < n; i++) {
251 fIndex[i] = t.fIndex[i];
252 fIndexBackup[i] = t.fIndex[i];
253 fdQdl[i] = t.fdQdl[i];
254 }
255 for (k = n; k < kMAXCLUSTERSPERTRACK; k++) {
256 fdQdl[k] = 0;
257 fIndex[k] = 0;
258 fIndexBackup[k] = 0;
259 }
260
261 for (i = 0; i < 6; i++) {
03b0452e 262 fTracklets[i] = t.fTracklets[i];
263 }
15ed8ba1 264
265 for (i = 0; i < 3; i++) {
266 fBudget[i] = t.fBudget[i];
267 }
268
5443e65e 269}
270
271//_____________________________________________________________________________
15ed8ba1 272AliTRDtrack::AliTRDtrack(const AliKalmanTrack &t, Double_t alpha)
273 :AliKalmanTrack(t)
274 ,fSeedLab(-1)
275 ,fdEdx(t.GetPIDsignal())
15ed8ba1 276 ,fDE(0)
277 ,fAlpha(alpha)
278 ,fX(0)
279 ,fStopped(kFALSE)
280 ,fY(0)
281 ,fZ(0)
282 ,fE(0)
283 ,fT(0)
284 ,fC(0)
285 ,fCyy(0)
286 ,fCzy(0)
287 ,fCzz(0)
288 ,fCey(0)
289 ,fCez(0)
290 ,fCee(0)
291 ,fCty(0)
292 ,fCtz(0)
293 ,fCte(0)
294 ,fCtt(0)
295 ,fCcy(0)
296 ,fCcz(0)
297 ,fCce(0)
298 ,fCct(0)
299 ,fCcc(0)
300 ,fLhElectron(0.0)
301 ,fNWrong(0)
302 ,fNRotate(0)
303 ,fNCross(0)
304 ,fNExpected(0)
305 ,fNLast(0)
306 ,fNExpectedLast(0)
307 ,fNdedx(0)
308 ,fChi2Last(0.0)
309 ,fBackupTrack(0x0)
53c17fbf 310{
5443e65e 311 //
312 // Constructor from AliTPCtrack or AliITStrack .
313 //
314
15ed8ba1 315 Int_t i = 0;
316 Int_t j = 0;
317 UInt_t k = 0;
318
319 SetChi2(0.0);
5443e65e 320 SetNumberOfClusters(0);
321
15ed8ba1 322 for (i = 0; i < kNplane; i++) {
323 for (j = 0; j < kNslice; j++) {
6d45eaef 324 fdEdxPlane[i][j] = 0.0;
325 }
eab5961e 326 fTimBinPlane[i] = -1;
327 }
5443e65e 328
15ed8ba1 329 if (fAlpha < -TMath::Pi()) {
330 fAlpha += 2.0 * TMath::Pi();
331 }
332 else if (fAlpha >= TMath::Pi()) {
333 fAlpha -= 2.0 * TMath::Pi();
334 }
79e94bf8 335
15ed8ba1 336 Double_t x;
337 Double_t p[5];
338 t.GetExternalParameters(x,p);
339 fX = x;
340 fY = p[0];
341 fZ = p[1];
342 fT = p[3];
343 x = GetLocalConvConst();
344 fC = p[4] / x;
345 fE = fC * fX - p[2];
346
347 // Conversion of the covariance matrix
348 Double_t c[15];
349 t.GetExternalCovariance(c);
350 c[10] /= x;
351 c[11] /= x;
352 c[12] /= x;
353 c[13] /= x;
354 c[14] /= x*x;
355
356 Double_t c22 = fX*fX * c[14] - 2.0*fX*c[12] + c[ 5];
357 Double_t c32 = fX * c[13] - c[ 8];
358 Double_t c20 = fX * c[10] - c[ 3];
359 Double_t c21 = fX * c[11] - c[ 4];
360 Double_t c42 = fX * c[14] - c[12];
361
362 fCyy = c[ 0];
363 fCzy = c[ 1]; fCzz = c[ 2];
364 fCey = c20; fCez = c21; fCee = c22;
365 fCty = c[ 6]; fCtz = c[ 7]; fCte = c32; fCtt = c[ 9];
366 fCcy = c[10]; fCcz = c[11]; fCce = c42; fCct = c[13]; fCcc = c[14];
367
368 for (k = 0; k < kMAXCLUSTERSPERTRACK; k++) {
369 fdQdl[k] = 0;
370 fIndex[k] = 0;
371 fIndexBackup[k] = 0;
79e94bf8 372 }
4f1c04d3 373
15ed8ba1 374 for (i = 0; i < 3; i++) {
375 fBudget[i] = 0;
376 }
377
79e94bf8 378}
53c17fbf 379
79e94bf8 380//_____________________________________________________________________________
15ed8ba1 381AliTRDtrack::AliTRDtrack(const AliESDtrack &t)
382 :AliKalmanTrack()
383 ,fSeedLab(-1)
384 ,fdEdx(t.GetTRDsignal())
15ed8ba1 385 ,fDE(0)
386 ,fAlpha(t.GetAlpha())
387 ,fX(0)
388 ,fStopped(kFALSE)
389 ,fY(0)
390 ,fZ(0)
391 ,fE(0)
392 ,fT(0)
393 ,fC(0)
394 ,fCyy(1e10)
395 ,fCzy(0)
396 ,fCzz(1e10)
397 ,fCey(0)
398 ,fCez(0)
399 ,fCee(1e10)
400 ,fCty(0)
401 ,fCtz(0)
402 ,fCte(0)
403 ,fCtt(1e10)
404 ,fCcy(0)
405 ,fCcz(0)
406 ,fCce(0)
407 ,fCct(0)
408 ,fCcc(1e10)
409 ,fLhElectron(0.0)
410 ,fNWrong(0)
411 ,fNRotate(0)
412 ,fNCross(0)
413 ,fNExpected(0)
414 ,fNLast(0)
415 ,fNExpectedLast(0)
416 ,fNdedx(0)
417 ,fChi2Last(0.0)
418 ,fBackupTrack(0x0)
53c17fbf 419{
79e94bf8 420 //
421 // Constructor from AliESDtrack
422 //
53c17fbf 423
15ed8ba1 424 Int_t i = 0;
425 Int_t j = 0;
426 UInt_t k = 0;
427
79e94bf8 428 SetLabel(t.GetLabel());
429 SetChi2(0.);
430 SetMass(t.GetMass());
1e9bb598 431 SetNumberOfClusters(t.GetTRDclusters(fIndex));
15ed8ba1 432
46e2d86c 433 Int_t ncl = t.GetTRDclusters(fIndexBackup);
15ed8ba1 434 for (k = ncl; k < kMAXCLUSTERSPERTRACK; k++) {
435 fIndexBackup[k] = 0;
436 fIndex[k] = 0;
437 }
438
439 for (i = 0; i < kNplane; i++) {
440 for (j = 0; j < kNslice; j++) {
6d45eaef 441 fdEdxPlane[i][j] = t.GetTRDsignals(i,j);
442 }
eab5961e 443 fTimBinPlane[i] = t.GetTRDTimBin(i);
444 }
79e94bf8 445
15ed8ba1 446 if (fAlpha < -TMath::Pi()) {
447 fAlpha += 2.0 * TMath::Pi();
448 }
449 else if (fAlpha >= TMath::Pi()) {
450 fAlpha -= 2.0 * TMath::Pi();
451 }
79e94bf8 452
15ed8ba1 453 // Conversion of the covariance matrix
454 Double_t x;
455 Double_t p[5];
456 t.GetExternalParameters(x,p);
457 Double_t c[15];
458 t.GetExternalCovariance(c);
459 if (t.GetStatus() & AliESDtrack::kTRDbackup) {
c0b978f0 460 t.GetOuterExternalParameters(fAlpha,x,p);
c9ec41e8 461 t.GetOuterExternalCovariance(c);
15ed8ba1 462 if (fAlpha < -TMath::Pi()) {
463 fAlpha += 2.0 * TMath::Pi();
464 }
465 else if (fAlpha >= TMath::Pi()) {
466 fAlpha -= 2.0 * TMath::Pi();
467 }
c5a8e3df 468 }
79e94bf8 469
15ed8ba1 470 fX = x;
471 fY = p[0];
472 fZ = p[1];
473 SaveLocalConvConst();
474 fT = p[3];
475 x = GetLocalConvConst();
476 fC = p[4] / x;
477 fE = fC*fX - p[2];
478
479 c[10] /= x;
480 c[11] /= x;
481 c[12] /= x;
482 c[13] /= x;
483 c[14] /= x*x;
484
485 Double_t c22 = fX*fX * c[14] - 2.0*fX*c[12] + c[ 5];
486 Double_t c32 = fX * c[13] - c[ 8];
487 Double_t c20 = fX * c[10] - c[ 3];
488 Double_t c21 = fX * c[11] - c[ 4];
489 Double_t c42 = fX * c[14] - c[12];
490
491 fCyy = c[ 0];
492 fCzy = c[ 1]; fCzz = c[ 2];
493 fCey = c20; fCez = c21; fCee = c22;
494 fCty = c[ 6]; fCtz = c[ 7]; fCte = c32; fCtt = c[ 9];
495 fCcy = c[10]; fCcz = c[11]; fCce = c42; fCct = c[13]; fCcc = c[14];
496
497 for (k = 0; k < kMAXCLUSTERSPERTRACK; k++) {
498 fdQdl[k] = 0;
499 //fIndex[k] = 0; //MI store indexes
500 }
5443e65e 501
15ed8ba1 502 for (i = 0; i < 3; i++) {
503 fBudget[i] = 0;
504 }
505 if ((t.GetStatus() & AliESDtrack::kTIME) == 0) {
506 return;
0d5b5c27 507 }
c630aafd 508
c630aafd 509 StartTimeIntegral();
15ed8ba1 510 Double_t times[10];
511 t.GetIntegratedTimes(times);
512 SetIntegratedTimes(times);
c630aafd 513 SetIntegratedLength(t.GetIntegratedLength());
514
16d9fbba 515}
516
53c17fbf 517//____________________________________________________________________________
518AliTRDtrack::~AliTRDtrack()
3fad3d32 519{
520 //
53c17fbf 521 // Destructor
522 //
3fad3d32 523
15ed8ba1 524 if (fBackupTrack) {
525 delete fBackupTrack;
526 }
53c17fbf 527 fBackupTrack = 0;
3fad3d32 528
53c17fbf 529}
530
531//____________________________________________________________________________
532AliTRDtrack &AliTRDtrack::operator=(const AliTRDtrack &t)
16d9fbba 533{
534 //
53c17fbf 535 // Assignment operator
16d9fbba 536 //
eab5961e 537
15ed8ba1 538 fLhElectron = 0.0;
539 fNWrong = 0;
540 fStopped = 0;
541 fNRotate = 0;
542 fNExpected = 0;
543 fNExpectedLast = 0;
544 fNdedx = 0;
545 fNCross = 0;
546 fNLast = 0;
547 fChi2Last = 0;
548 fBackupTrack = 0;
53c17fbf 549
550 fAlpha = t.GetAlpha();
15ed8ba1 551 if (fAlpha < -TMath::Pi()) {
552 fAlpha += 2.0 * TMath::Pi();
553 }
554 else if (fAlpha >= TMath::Pi()) {
555 fAlpha -= 2.0 * TMath::Pi();
556 }
53c17fbf 557
558 return *this;
eab5961e 559
16d9fbba 560}
9c9d2487 561
53c17fbf 562//____________________________________________________________________________
563Float_t AliTRDtrack::StatusForTOF()
7ad19338 564{
53c17fbf 565 //
566 // Defines the status of the TOF extrapolation
567 //
03b0452e 568
15ed8ba1 569 // Definition of res ????
570 Float_t res = (0.2 + 0.8 * (fN / (fNExpected + 5.0)))
571 * (0.4 + 0.6 * fTracklets[5].GetN() / 20.0);
572 res *= (0.25 + 0.8 * 40.0 / (40.0 + fBudget[2]));
03b0452e 573 return res;
574
15ed8ba1 575 // This part of the function is never reached ????
576 // What defines these parameters ????
577 Int_t status = 0;
578 if (GetNumberOfClusters() < 20) return 0;
579 if ((fN > 110) &&
580 (fChi2/(Float_t(fN)) < 3.0)) return 3; // Gold
581 if ((fNLast > 30) &&
582 (fChi2Last/(Float_t(fNLast)) < 3.0)) return 3; // Gold
583 if ((fNLast > 20) &&
584 (fChi2Last/(Float_t(fNLast)) < 2.0)) return 3; // Gold
585 if ((fNLast/(fNExpectedLast+3.0) > 0.8) &&
586 (fChi2Last/Float_t(fNLast) < 5.0) &&
587 (fNLast > 20)) return 2; // Silver
588 if ((fNLast > 5) &&
589 (((fNLast+1.0)/(fNExpectedLast+1.0)) > 0.8) &&
590 (fChi2Last/(fNLast-5.0) < 6.0)) return 1;
53c17fbf 591
4f1c04d3 592 return status;
53c17fbf 593
4f1c04d3 594}
595
5443e65e 596//_____________________________________________________________________________
53c17fbf 597void AliTRDtrack::GetExternalCovariance(Double_t cc[15]) const
598{
5443e65e 599 //
600 // This function returns external representation of the covriance matrix.
601 //
53c17fbf 602
15ed8ba1 603 Double_t a = GetLocalConvConst();
5443e65e 604
15ed8ba1 605 Double_t c22 = fX*fX*fCcc - 2.0*fX*fCce+fCee;
606 Double_t c32 = fX*fCct-fCte;
607 Double_t c20 = fX*fCcy-fCey;
608 Double_t c21 = fX*fCcz-fCez;
609 Double_t c42 = fX*fCcc-fCce;
5443e65e 610
15ed8ba1 611 cc[ 0]=fCyy;
612 cc[ 1]=fCzy; cc[ 2]=fCzz;
613 cc[ 3]=c20; cc[ 4]=c21; cc[ 5]=c22;
614 cc[ 6]=fCty; cc[ 7]=fCtz; cc[ 8]=c32; cc[ 9]=fCtt;
b8dc2353 615 cc[10]=fCcy*a; cc[11]=fCcz*a; cc[12]=c42*a; cc[13]=fCct*a; cc[14]=fCcc*a*a;
616
5443e65e 617}
618
46d29e70 619//_____________________________________________________________________________
53c17fbf 620void AliTRDtrack::GetCovariance(Double_t cc[15]) const
621{
622 //
623 // Returns the track covariance matrix
624 //
b8dc2353 625
15ed8ba1 626 cc[ 0]=fCyy;
627 cc[ 1]=fCzy; cc[ 2]=fCzz;
628 cc[ 3]=fCey; cc[ 4]=fCez; cc[ 5]=fCee;
629 cc[ 6]=fCcy; cc[ 7]=fCcz; cc[ 8]=fCce; cc[ 9]=fCcc;
b3a5a838 630 cc[10]=fCty; cc[11]=fCtz; cc[12]=fCte; cc[13]=fCct; cc[14]=fCtt;
b8dc2353 631
46d29e70 632}
633
634//_____________________________________________________________________________
53c17fbf 635Int_t AliTRDtrack::Compare(const TObject *o) const
636{
637 //
638 // Compares tracks according to their Y2 or curvature
639 //
46d29e70 640
d49e463b 641 AliTRDtrack *t = (AliTRDtrack *) o;
46d29e70 642 // Double_t co=t->GetSigmaY2();
643 // Double_t c =GetSigmaY2();
644
d49e463b 645 Double_t co = TMath::Abs(t->GetC());
646 Double_t c = TMath::Abs(GetC());
46d29e70 647
d49e463b 648 if (c > co) {
649 return 1;
650 }
651 else if (c < co) {
652 return -1;
653 }
46d29e70 654 return 0;
53c17fbf 655
46d29e70 656}
657
a819a5f7 658//_____________________________________________________________________________
15ed8ba1 659void AliTRDtrack::CookdEdx(Double_t low, Double_t up)
660{
661 //
d49e463b 662 // Calculates the truncated dE/dx within the "low" and "up" cuts.
15ed8ba1 663 //
664
d49e463b 665 Int_t i = 0;
666
667 // Array to sort the dEdx values according to amplitude
668 Float_t sorted[kMAXCLUSTERSPERTRACK];
15ed8ba1 669
670 // Number of clusters used for dedx
d49e463b 671 Int_t nc = fNdedx;
672
15ed8ba1 673 // Require at least 10 clusters for a dedx measurement
15ed8ba1 674 if (nc < 10) {
4f1c04d3 675 SetdEdx(0);
676 return;
677 }
a819a5f7 678
15ed8ba1 679 // Lower and upper bound
680 Int_t nl = Int_t(low * nc);
681 Int_t nu = Int_t( up * nc);
682
d49e463b 683 // Can fdQdl be negative ????
15ed8ba1 684 for (i = 0; i < nc; i++) {
685 sorted[i] = TMath::Abs(fdQdl[i]);
686 }
687
688 // Sort the dedx values by amplitude
689 Int_t *index = new Int_t[nc];
690 TMath::Sort(nc,sorted,index,kFALSE);
d49e463b 691
15ed8ba1 692 // Sum up the truncated charge between nl and nu
d49e463b 693 Float_t dedx = 0.0;
15ed8ba1 694 for (i = nl; i <= nu; i++) {
695 dedx += sorted[index[i]];
696 }
d49e463b 697 dedx /= (nu - nl + 1.0);
3551db50 698 SetdEdx(dedx);
699
15ed8ba1 700 delete [] index;
701
a819a5f7 702}
703
46d29e70 704//_____________________________________________________________________________
b3a5a838 705Int_t AliTRDtrack::PropagateTo(Double_t xk,Double_t x0,Double_t rho)
46d29e70 706{
15ed8ba1 707 //
46d29e70 708 // Propagates a track of particle with mass=pm to a reference plane
709 // defined by x=xk through media of density=rho and radiationLength=x0
15ed8ba1 710 //
711
712 if (xk == fX) {
713 return 1;
714 }
715
716 if (TMath::Abs(fC*xk - fE) >= 0.9) {
717 return 0;
718 }
719
720 Double_t lcc = GetLocalConvConst();
46d29e70 721
15ed8ba1 722 Double_t oldX = fX;
723 Double_t oldY = fY;
724 Double_t oldZ = fZ;
9c9d2487 725
15ed8ba1 726 Double_t x1 = fX;
727 Double_t x2 = x1 + (xk - x1);
728 Double_t dx = x2 - x1;
729 Double_t y1 = fY;
730 Double_t z1 = fZ;
731 Double_t c1 = fC*x1 - fE;
732 if ((c1*c1) > 1) {
46d29e70 733 return 0;
734 }
15ed8ba1 735
736 Double_t r1 = TMath::Sqrt(1.0 - c1*c1);
737 Double_t c2 = fC*x2 - fE;
738 if ((c2*c2) > 1) {
739 return 0;
740 }
741
742 Double_t r2 = TMath::Sqrt(1.0 - c2*c2);
743
744 fY += dx*(c1+c2) / (r1+r2);
745 fZ += dx*(c1+c2) / (c1*r2 + c2*r1) * fT;
746
747 // f = F - 1
748 Double_t rr = r1+r2;
749 Double_t cc = c1+c2;
750 Double_t xx = x1+x2;
751 Double_t f02 = -dx*(2*rr + cc*(c1/r1 + c2/r2))/(rr*rr);
752 Double_t f04 = dx*(rr*xx + cc*(c1*x1/r1+c2*x2/r2))/(rr*rr);
753 Double_t cr = c1*r2+c2*r1;
754 Double_t f12 = -dx*fT*(2*cr + cc*(c2*c1/r1-r1 + c1*c2/r2-r2))/(cr*cr);
755 Double_t f13 = dx*cc/cr;
756 Double_t f14 = dx*fT*(cr*xx-cc*(r1*x2-c2*c1*x1/r1+r2*x1-c1*c2*x2/r2))/(cr*cr);
757
758 // b = C*ft
759 Double_t b00 = f02*fCey + f04*fCcy;
760 Double_t b01 = f12*fCey + f14*fCcy + f13*fCty;
761 Double_t b10 = f02*fCez + f04*fCcz;
762 Double_t b11 = f12*fCez + f14*fCcz + f13*fCtz;
763 Double_t b20 = f02*fCee + f04*fCce;
764 Double_t b21 = f12*fCee + f14*fCce + f13*fCte;
765 Double_t b30 = f02*fCte + f04*fCct;
766 Double_t b31 = f12*fCte + f14*fCct + f13*fCtt;
767 Double_t b40 = f02*fCce + f04*fCcc;
768 Double_t b41 = f12*fCce + f14*fCcc + f13*fCct;
769
770 // a = f*b = f*C*ft
771 Double_t a00 = f02*b20 + f04*b40;
772 Double_t a01 = f02*b21 + f04*b41;
773 Double_t a11 = f12*b21 + f14*b41 + f13*b31;
774
775 // F*C*Ft = C + (a + b + bt)
776 fCyy += a00 + 2.0*b00;
46d29e70 777 fCzy += a01 + b01 + b10;
b3a5a838 778 fCey += b20;
779 fCty += b30;
780 fCcy += b40;
15ed8ba1 781 fCzz += a11 + 2.0*b11;
b3a5a838 782 fCez += b21;
783 fCtz += b31;
b8dc2353 784 fCcz += b41;
46d29e70 785
15ed8ba1 786 fX = x2;
46d29e70 787
15ed8ba1 788 // Change of the magnetic field
c84a5e9e 789 SaveLocalConvConst();
15ed8ba1 790 cc = fC;
791 fC *= lcc / GetLocalConvConst();
792 fE += fX * (fC-cc);
793
794 // Multiple scattering
795 // What is 14.1 ????
796 Double_t d = TMath::Sqrt((x1-fX)*(x1-fX) + (y1-fY)*(y1-fY) + (z1-fZ)*(z1-fZ));
797 Double_t p2 = (1.0 + GetTgl()*GetTgl()) / (Get1Pt()*Get1Pt());
798 Double_t beta2 = p2 / (p2 + GetMass()*GetMass());
799 Double_t theta2 = 14.1*14.1 / (beta2*p2*1e6) * d / x0 * rho;
800 Double_t ey = fC*fX - fE;
801 Double_t ez = fT;
802 Double_t xz = fC*ez;
803 Double_t zz1 = ez*ez + 1.0;
804 Double_t xy = fE + ey;
3c625a9b 805
15ed8ba1 806 fCee += (2.0*ey*ez*ez*fE + 1.0 - ey*ey + ez*ez + fE*fE*ez*ez) * theta2;
46d29e70 807 fCte += ez*zz1*xy*theta2;
808 fCtt += zz1*zz1*theta2;
b3a5a838 809 fCce += xz*ez*xy*theta2;
810 fCct += xz*zz1*theta2;
811 fCcc += xz*xz*theta2;
5443e65e 812
15ed8ba1 813 // Energy losses
814 // What is 5940.0 ???? and 0.153e-3 ????
815 if ((5940.0*beta2 / (1.0 - beta2 + 1e-10) - beta2) < 0.0) {
816 return 0;
817 }
818 Double_t dE = 0.153e-3/beta2 * (TMath::Log(5940.0*beta2 / (1.0 - beta2 + 1e-10)) - beta2)
819 * d * rho;
820 Float_t budget = d * rho;
2276b464 821 fBudget[0] +=budget;
15ed8ba1 822
823 // Suspicious part - think about it ????
96c3a73c 824 Double_t kinE = TMath::Sqrt(p2);
15ed8ba1 825 if (dE > 0.8 * kinE) {
826 dE = 0.8 * kinE;
827 }
828 if (dE < 0.0) {
829 dE = 0.0; // Not valid region for Bethe bloch
830 }
831 fDE += dE;
832 if (x1 < x2) {
833 dE = -dE;
834 }
835 cc = fC;
836 fC *= (1.0 - TMath::Sqrt(p2 + GetMass()*GetMass()) / p2 * dE);
837 fE += fX * (fC - cc);
838
839 // Energy loss fluctuation
840 // Why 0.07 ????
841 Double_t sigmade = 0.07 * TMath::Sqrt(TMath::Abs(dE));
842 Double_t sigmac2 = sigmade*sigmade * fC*fC * (p2 + GetMass()*GetMass()) / (p2*p2);
3fad3d32 843 fCcc += sigmac2;
15ed8ba1 844 fCee += fX*fX * sigmac2;
845
846 // Track time measurement
847 if (x1 < x2) {
848 if (IsStartedTimeIntegral()) {
849 Double_t l2 = TMath::Sqrt((fX-oldX)*(fX-oldX)
850 + (fY-oldY)*(fY-oldY)
851 + (fZ-oldZ)*(fZ-oldZ));
852 if (TMath::Abs(l2*fC) > 0.0001){
853 // <ake correction for curvature if neccesary
854 l2 = 0.5 * TMath::Sqrt((fX-oldX)*(fX-oldX)
855 + (fY-oldY)*(fY-oldY));
856 l2 = 2.0 * TMath::ASin(l2 * fC) / fC;
857 l2 = TMath::Sqrt(l2*l2 + (fZ-oldZ)*(fZ-oldZ));
858 }
859 AddTimeStep(l2);
03b0452e 860 }
0d5b5c27 861 }
862
b8dc2353 863 return 1;
6d45eaef 864
46d29e70 865}
866
46d29e70 867//_____________________________________________________________________________
15ed8ba1 868Int_t AliTRDtrack::Update(const AliTRDcluster *c, Double_t chisq
869 , UInt_t index, Double_t h01)
46d29e70 870{
15ed8ba1 871 //
872 // Assignes a found cluster to the track and updates track information
873 //
46d29e70 874
b8dc2353 875 Bool_t fNoTilt = kTRUE;
15ed8ba1 876 // What is 0.003 ????
877 if (TMath::Abs(h01) > 0.003) {
878 fNoTilt = kFALSE;
46e2d86c 879 }
46d29e70 880
15ed8ba1 881 // Add angular effect to the error contribution
882 Float_t tangent2 = (fC*fX-fE) * (fC*fX-fE);
883 if (tangent2 < 0.90000){
884 tangent2 = tangent2 / (1.0 - tangent2);
885 }
886 // What is 0.04 ????
887 Float_t errang = tangent2 * 0.04;
888 Float_t padlength = TMath::Sqrt(c->GetSigmaZ2() * 12.0);
889
890 Double_t r00 = c->GetSigmaY2() + errang;
891 Double_t r01 = 0.0;
892 Double_t r11 = c->GetSigmaZ2() * 100.0;
893 r00 += fCyy;
894 r01 += fCzy;
895 r11 += fCzz;
896 Double_t det = r00*r11 - r01*r01;
897 Double_t tmp = r00;
898 r00 = r11 / det;
899 r11 = tmp / det;
900 r01 = -r01 / det;
901
902 Double_t k00 = fCyy*r00 + fCzy*r01;
903 Double_t k01 = fCyy*r01 + fCzy*r11;
904 Double_t k10 = fCzy*r00 + fCzz*r01;
905 Double_t k11 = fCzy*r01 + fCzz*r11;
906 Double_t k20 = fCey*r00 + fCez*r01;
907 Double_t k21 = fCey*r01 + fCez*r11;
908 Double_t k30 = fCty*r00 + fCtz*r01;
909 Double_t k31 = fCty*r01 + fCtz*r11;
910 Double_t k40 = fCcy*r00 + fCcz*r01;
911 Double_t k41 = fCcy*r01 + fCcz*r11;
912
913 Double_t dy = c->GetY() - fY;
914 Double_t dz = c->GetZ() - fZ;
915 Double_t cur = fC + k40*dy + k41*dz;
916 Double_t eta = fE + k20*dy + k21*dz;
917
918 if (fNoTilt) {
919
920 if (TMath::Abs(cur*fX-eta) >= 0.9) {
b8dc2353 921 return 0;
922 }
923 fY += k00*dy + k01*dz;
924 fZ += k10*dy + k11*dz;
925 fE = eta;
b8dc2353 926 fC = cur;
15ed8ba1 927
b8dc2353 928 }
929 else {
b8dc2353 930
15ed8ba1 931 // Empirical factor set by C.Xu in the first tilt version
932 // Is this factor still ok ????
933 Double_t xuFactor = 100.0;
934 dy = c->GetY() - fY;
935 dz = c->GetZ() - fZ;
936 dy = dy + h01*dz;
937
938 Float_t add = 0.0;
939 if (TMath::Abs(dz) > padlength/2.0) {
940 Float_t dy2 = c->GetY() - fY;
941 Float_t sign = (dz > 0.0) ? -1.0 : 1.0;
942 dy2 += h01 * sign * padlength/2.0;
943 dy = dy2;
944 add = 0.0;
945 }
946
947 r00 = c->GetSigmaY2() + errang + add;
948 r01 = 0.0;
949 r11 = c->GetSigmaZ2() * xuFactor;
950 r00 += (fCyy + 2.0*h01*fCzy + h01*h01*fCzz);
951 r01 += (fCzy + h01*fCzz);
952 r11 += fCzz;
953
954 det = r00*r11 - r01*r01;
955 tmp = r00;
956 r00 = r11/det;
957 r11 = tmp/det;
958 r01 = -r01/det;
959
960 k00 = fCyy*r00 + fCzy * (r01 + h01*r00);
961 k01 = fCyy*r01 + fCzy * (r11 + h01*r01);
962 k10 = fCzy*r00 + fCzz * (r01 + h01*r00);
963 k11 = fCzy*r01 + fCzz * (r11 + h01*r01);
964 k20 = fCey*r00 + fCez * (r01 + h01*r00);
965 k21 = fCey*r01 + fCez * (r11 + h01*r01);
966 k30 = fCty*r00 + fCtz * (r01 + h01*r00);
967 k31 = fCty*r01 + fCtz * (r11 + h01*r01);
968 k40 = fCcy*r00 + fCcz * (r01 + h01*r00);
969 k41 = fCcy*r01 + fCcz * (r11 + h01*r01);
970
971 cur = fC + k40*dy + k41*dz;
972 eta = fE + k20*dy + k21*dz;
973 if (TMath::Abs(cur*fX - eta) >= 0.9) {
b8dc2353 974 return 0;
975 }
15ed8ba1 976
977 fY += k00*dy + k01*dz;
978 fZ += k10*dy + k11*dz;
979 fE = eta;
980 fT += k30*dy + k31*dz;
981 fC = cur;
b8dc2353 982
15ed8ba1 983 k01 += h01*k00;
984 k11 += h01*k10;
985 k21 += h01*k20;
986 k31 += h01*k30;
987 k41 += h01*k40;
46e2d86c 988
b8dc2353 989 }
46e2d86c 990
15ed8ba1 991 Double_t c01 = fCzy;
992 Double_t c02 = fCey;
993 Double_t c03 = fCty;
994 Double_t c04 = fCcy;
995 Double_t c12 = fCez;
996 Double_t c13 = fCtz;
997 Double_t c14 = fCcz;
998
999 fCyy -= k00*fCyy + k01*fCzy;
1000 fCzy -= k00*c01 + k01*fCzz;
1001 fCey -= k00*c02 + k01*c12;
1002 fCty -= k00*c03 + k01*c13;
1003 fCcy -= k00*c04 + k01*c14;
b8dc2353 1004
15ed8ba1 1005 fCzz -= k10*c01 + k11*fCzz;
1006 fCez -= k10*c02 + k11*c12;
1007 fCtz -= k10*c03 + k11*c13;
1008 fCcz -= k10*c04 + k11*c14;
b8dc2353 1009
15ed8ba1 1010 fCee -= k20*c02 + k21*c12;
1011 fCte -= k20*c03 + k21*c13;
1012 fCce -= k20*c04 + k21*c14;
b8dc2353 1013
15ed8ba1 1014 fCtt -= k30*c03 + k31*c13;
1015 fCct -= k40*c03 + k41*c13;
b8dc2353 1016
15ed8ba1 1017 fCcc -= k40*c04 + k41*c14;
46d29e70 1018
15ed8ba1 1019 Int_t n = GetNumberOfClusters();
1020 fIndex[n] = index;
b8dc2353 1021 SetNumberOfClusters(n+1);
fd621f36 1022
b8dc2353 1023 SetChi2(GetChi2()+chisq);
5443e65e 1024
b8dc2353 1025 return 1;
53c17fbf 1026
46d29e70 1027}
53c17fbf 1028
46e2d86c 1029//_____________________________________________________________________________
15ed8ba1 1030Int_t AliTRDtrack::UpdateMI(const AliTRDcluster *c, Double_t chisq
1031 , UInt_t index, Double_t h01
1032 , Int_t /*plane*/)
46e2d86c 1033{
15ed8ba1 1034 //
46e2d86c 1035 // Assignes found cluster to the track and updates track information
15ed8ba1 1036 //
46e2d86c 1037
1038 Bool_t fNoTilt = kTRUE;
15ed8ba1 1039 if (TMath::Abs(h01) > 0.003) {
1040 fNoTilt = kFALSE;
1041 }
1042
1043 //
1044 // Add angular effect to the error contribution and make correction
1045 // Still needed ????
1046 // AliTRDclusterCorrection *corrector = AliTRDclusterCorrection::GetCorrection();
3c625a9b 1047 //
4f1c04d3 1048
15ed8ba1 1049 Double_t tangent2 = (fC*fX-fE) * (fC*fX-fE);
1050 if (tangent2 < 0.9) {
1051 tangent2 = tangent2 / (1.0 - tangent2);
1052 }
1053 Double_t tangent = TMath::Sqrt(tangent2);
1054 if ((fC*fX-fE) < 0.0) {
1055 tangent *= -1.0;
1056 }
1057
1058 // Where are the parameters from ????
1059 Double_t errang = tangent2 * 0.04;
1060 Double_t errsys = 0.025*0.025 * 20.0; // Systematic error part
1061
1062 Float_t extend = 1.0;
1063 if (c->GetNPads() == 4) extend = 2.0;
1064
1065 /////////////////////////////////////////////////////////////////////////////
1066 //
1067 // Is this still needed or will it be needed ????
1068 //
1069 //if (c->GetNPads() == 5) extend = 3.0;
1070 //if (c->GetNPads() == 6) extend = 3.0;
1071 //if (c->GetQ() < 15) {
1072 // return 1;
1073 //}
1074 //
1075 // Will this be needed ????
c5a8e3df 1076 /*
15ed8ba1 1077 if (corrector !=0 ) {
3c625a9b 1078 //if (0){
1079 correction = corrector->GetCorrection(plane,c->GetLocalTimeBin(),tangent);
1080 if (TMath::Abs(correction)>0){
1081 //if we have info
1082 errang = corrector->GetSigma(plane,c->GetLocalTimeBin(),tangent);
1083 errang *= errang;
1084 errang += tangent2*0.04;
1085 }
1086 }
c5a8e3df 1087 */
15ed8ba1 1088 //
c84a5e9e 1089 // Double_t padlength = TMath::Sqrt(c->GetSigmaZ2()*12.);
15ed8ba1 1090 /////////////////////////////////////////////////////////////////////////////
1091
1092 Double_t r00 = (c->GetSigmaY2() + errang + errsys) * extend;
1093 Double_t r01 = 0.0;
1094 Double_t r11 = c->GetSigmaZ2()*10000.0;
1095 r00 += fCyy;
1096 r01 += fCzy;
1097 r11 += fCzz;
1098 Double_t det =r00*r11 - r01*r01;
1099 Double_t tmp =r00;
1100 r00 = r11 / det;
1101 r11 = tmp / det;
1102 r01 = -r01 / det;
1103
1104 Double_t k00 = fCyy*r00 + fCzy*r01;
1105 Double_t k01 = fCyy*r01 + fCzy*r11;
1106 Double_t k10 = fCzy*r00 + fCzz*r01;
1107 Double_t k11 = fCzy*r01 + fCzz*r11;
1108 Double_t k20 = fCey*r00 + fCez*r01;
1109 Double_t k21 = fCey*r01 + fCez*r11;
1110 Double_t k30 = fCty*r00 + fCtz*r01;
1111 Double_t k31 = fCty*r01 + fCtz*r11;
1112 Double_t k40 = fCcy*r00 + fCcz*r01;
1113 Double_t k41 = fCcy*r01 + fCcz*r11;
1114
1115 Double_t dy = c->GetY() - fY;
1116 Double_t dz = c->GetZ() - fZ;
1117 Double_t cur = fC + k40*dy + k41*dz;
1118 Double_t eta = fE + k20*dy + k21*dz;
1119
1120 if (fNoTilt) {
1121
1122 if (TMath::Abs(cur*fX - eta) >= 0.9) {
46e2d86c 1123 return 0;
1124 }
15ed8ba1 1125
46e2d86c 1126 fY += k00*dy + k01*dz;
1127 fZ += k10*dy + k11*dz;
1128 fE = eta;
46e2d86c 1129 fC = cur;
15ed8ba1 1130
46e2d86c 1131 }
1132 else {
15ed8ba1 1133
1134 Double_t padlength = TMath::Sqrt(c->GetSigmaZ2() * 12.0);
1135 // Empirical factor set by C.Xu in the first tilt version
1136 Double_t xuFactor = 1000.0;
1137
1138 dy = c->GetY() - fY;
1139 dz = c->GetZ() - fZ;
1140 //dy = dy + h01*dz + correction; // Still needed ????
4f1c04d3 1141
1142 Double_t tiltdz = dz;
15ed8ba1 1143 if (TMath::Abs(tiltdz) > padlength/2.0) {
1144 tiltdz = TMath::Sign(padlength/2.0,dz);
c84a5e9e 1145 }
15ed8ba1 1146 dy = dy + h01*tiltdz;
4f1c04d3 1147
15ed8ba1 1148 Double_t add = 0.0;
1149 if (TMath::Abs(dz) > padlength/2.0) {
1150 //Double_t dy2 = c->GetY() - fY; // Still needed ????
c84a5e9e 1151 //Double_t sign = (dz>0) ? -1.: 1.;
1152 //dy2-=h01*sign*padlength/2.;
1153 //dy = dy2;
15ed8ba1 1154 add = 1.0;
c84a5e9e 1155 }
15ed8ba1 1156
1157 Double_t s00 = (c->GetSigmaY2() + errang) * extend + errsys + add; // Error pad
1158 Double_t s11 = c->GetSigmaZ2() * xuFactor; // Error pad-row
46e2d86c 1159 //
15ed8ba1 1160 r00 = fCyy + 2*fCzy*h01 + fCzz*h01*h01 + s00;
1161 r01 = fCzy + fCzz*h01;
1162 r11 = fCzz + s11;
1163 det = r00*r11 - r01*r01;
1164
1165 // Inverse matrix
1166 tmp = r00;
1167 r00 = r11 / det;
1168 r11 = tmp / det;
1169 r01 = -r01 / det;
46e2d86c 1170
1171 // K matrix
15ed8ba1 1172 k00 = fCyy*r00 + fCzy * (r01 + h01*r00);
1173 k01 = fCyy*r01 + fCzy * (r11 + h01*r01);
1174 k10 = fCzy*r00 + fCzz * (r01 + h01*r00);
1175 k11 = fCzy*r01 + fCzz * (r11 + h01*r01);
1176 k20 = fCey*r00 + fCez * (r01 + h01*r00);
1177 k21 = fCey*r01 + fCez * (r11 + h01*r01);
1178 k30 = fCty*r00 + fCtz * (r01 + h01*r00);
1179 k31 = fCty*r01 + fCtz * (r11 + h01*r01);
1180 k40 = fCcy*r00 + fCcz * (r01 + h01*r00);
1181 k41 = fCcy*r01 + fCcz * (r11 + h01*r01);
1182
1183 // Update measurement
1184 cur = fC + k40*dy + k41*dz;
1185 eta = fE + k20*dy + k21*dz;
1186 if (TMath::Abs(cur*fX - eta) >= 0.9) {
46e2d86c 1187 return 0;
1188 }
15ed8ba1 1189 fY += k00*dy + k01*dz;
1190 fZ += k10*dy + k11*dz;
1191 fE = eta;
1192 fT += k30*dy + k31*dz;
1193 fC = cur;
46e2d86c 1194
15ed8ba1 1195 k01 += h01*k00;
1196 k11 += h01*k10;
1197 k21 += h01*k20;
1198 k31 += h01*k30;
1199 k41 += h01*k40;
46e2d86c 1200
1201 }
15ed8ba1 1202
1203 // Update the covariance matrix
1204 Double_t oldyy = fCyy;
1205 Double_t oldzz = fCzz;
1206 //Double_t oldee = fCee;
1207 //Double_t oldcc = fCcc;
1208 Double_t oldzy = fCzy;
1209 Double_t oldey = fCey;
1210 Double_t oldty = fCty;
1211 Double_t oldcy = fCcy;
1212 Double_t oldez = fCez;
1213 Double_t oldtz = fCtz;
1214 Double_t oldcz = fCcz;
1215 //Double_t oldte = fCte;
1216 //Double_t oldce = fCce;
46e2d86c 1217 //Double_t oldct = fCct;
1218
15ed8ba1 1219 fCyy -= k00*oldyy + k01*oldzy;
1220 fCzy -= k10*oldyy + k11*oldzy;
1221 fCey -= k20*oldyy + k21*oldzy;
1222 fCty -= k30*oldyy + k31*oldzy;
1223 fCcy -= k40*oldyy + k41*oldzy;
1224
1225 fCzz -= k10*oldzy + k11*oldzz;
1226 fCez -= k20*oldzy + k21*oldzz;
1227 fCtz -= k30*oldzy + k31*oldzz;
1228 fCcz -= k40*oldzy + k41*oldzz;
1229
1230 fCee -= k20*oldey + k21*oldez;
1231 fCte -= k30*oldey + k31*oldez;
1232 fCce -= k40*oldey + k41*oldez;
1233
1234 fCtt -= k30*oldty + k31*oldtz;
1235 fCct -= k40*oldty + k41*oldtz;
1236
1237 fCcc -= k40*oldcy + k41*oldcz;
46e2d86c 1238
15ed8ba1 1239 Int_t n = GetNumberOfClusters();
1240 fIndex[n] = index;
46e2d86c 1241 SetNumberOfClusters(n+1);
1242
15ed8ba1 1243 SetChi2(GetChi2() + chisq);
46e2d86c 1244
1245 return 1;
7ad19338 1246
53c17fbf 1247}
7ad19338 1248
7ad19338 1249//_____________________________________________________________________________
1250Int_t AliTRDtrack::UpdateMI(const AliTRDtracklet &tracklet)
1251{
1252 //
1253 // Assignes found tracklet to the track and updates track information
1254 //
15ed8ba1 1255
1256 Double_t r00 = (tracklet.GetTrackletSigma2());
1257 Double_t r01 = 0.0;
1258 Double_t r11 = 10000.0;
1259 r00 += fCyy;
1260 r01 += fCzy;
1261 r11 += fCzz;
1262
1263 Double_t det = r00*r11 - r01*r01;
1264 Double_t tmp = r00;
1265 r00 = r11 / det;
1266 r11 = tmp / det;
1267 r01 = -r01 / det;
7ad19338 1268
15ed8ba1 1269 Double_t dy = tracklet.GetY() - fY;
1270 Double_t dz = tracklet.GetZ() - fZ;
7ad19338 1271
15ed8ba1 1272 Double_t s00 = tracklet.GetTrackletSigma2(); // Error pad
1273 Double_t s11 = 100000.0; // Error pad-row
7ad19338 1274 Float_t h01 = tracklet.GetTilt();
15ed8ba1 1275
1276 r00 = fCyy + fCzz*h01*h01 + s00;
1277 r01 = fCzy;
7ad19338 1278 r11 = fCzz + s11;
1279 det = r00*r11 - r01*r01;
15ed8ba1 1280
1281 // Inverse matrix
1282 tmp = r00;
1283 r00 = r11 / det;
1284 r11 = tmp / det;
1285 r01 = -r01 / det;
1286
7ad19338 1287 // K matrix
15ed8ba1 1288 Double_t k00 = fCyy*r00 + fCzy*r01;
1289 Double_t k01 = fCyy*r01 + fCzy*r11;
1290 Double_t k10 = fCzy*r00 + fCzz*r01;
1291 Double_t k11 = fCzy*r01 + fCzz*r11;
1292 Double_t k20 = fCey*r00 + fCez*r01;
1293 Double_t k21 = fCey*r01 + fCez*r11;
1294 Double_t k30 = fCty*r00 + fCtz*r01;
1295 Double_t k31 = fCty*r01 + fCtz*r11;
1296 Double_t k40 = fCcy*r00 + fCcz*r01;
1297 Double_t k41 = fCcy*r01 + fCcz*r11;
1298
1299 // Update measurement
1300 Double_t cur = fC + k40*dy + k41*dz;
1301 Double_t eta = fE + k20*dy + k21*dz;
1302 if (TMath::Abs(cur*fX-eta) >= 0.9) {
7ad19338 1303 return 0;
1304 }
7ad19338 1305
1306 fY += k00*dy + k01*dz;
1307 fZ += k10*dy + k11*dz;
1308 fE = eta;
1309 fT += k30*dy + k31*dz;
1310 fC = cur;
1311
15ed8ba1 1312 // Update the covariance matrix
1313 Double_t oldyy = fCyy;
1314 Double_t oldzz = fCzz;
1315 //Double_t oldee = fCee;
1316 //Double_t oldcc = fCcc;
1317 Double_t oldzy = fCzy;
1318 Double_t oldey = fCey;
1319 Double_t oldty = fCty;
1320 Double_t oldcy = fCcy;
1321 Double_t oldez = fCez;
1322 Double_t oldtz = fCtz;
1323 Double_t oldcz = fCcz;
1324 //Double_t oldte = fCte;
1325 //Double_t oldce = fCce;
7ad19338 1326 //Double_t oldct = fCct;
1327
15ed8ba1 1328 fCyy -= k00*oldyy + k01*oldzy;
1329 fCzy -= k10*oldyy + k11*oldzy;
1330 fCey -= k20*oldyy + k21*oldzy;
1331 fCty -= k30*oldyy + k31*oldzy;
1332 fCcy -= k40*oldyy + k41*oldzy;
1333
1334 fCzz -= k10*oldzy + k11*oldzz;
1335 fCez -= k20*oldzy + k21*oldzz;
1336 fCtz -= k30*oldzy + k31*oldzz;
1337 fCcz -= k40*oldzy + k41*oldzz;
1338
1339 fCee -= k20*oldey + k21*oldez;
1340 fCte -= k30*oldey + k31*oldez;
1341 fCce -= k40*oldey + k41*oldez;
1342
1343 fCtt -= k30*oldty + k31*oldtz;
1344 fCct -= k40*oldty + k41*oldtz;
1345
1346 fCcc -= k40*oldcy + k41*oldcz;
7ad19338 1347
53c17fbf 1348 return 1;
7ad19338 1349
53c17fbf 1350}
c84a5e9e 1351
46d29e70 1352//_____________________________________________________________________________
3fad3d32 1353Int_t AliTRDtrack::Rotate(Double_t alpha, Bool_t absolute)
46d29e70 1354{
15ed8ba1 1355 //
1356 // Rotates the track parameters in the R*phi plane,
1357 // if the absolute rotation alpha is in global system.
1358 // Otherwise the rotation is relative to the current rotation angle
1359 //
1360
3fad3d32 1361 if (absolute) {
1362 alpha -= fAlpha;
1363 }
1364 else{
1365 fNRotate++;
1366 }
46d29e70 1367
1368 fAlpha += alpha;
15ed8ba1 1369 if (fAlpha < -TMath::Pi()) {
1370 fAlpha += 2.0 * TMath::Pi();
46d29e70 1371 }
15ed8ba1 1372 if (fAlpha >= TMath::Pi()) {
1373 fAlpha -= 2.0 * TMath::Pi();
1374 }
1375
1376 Double_t x1 = fX;
1377 Double_t y1 = fY;
1378 Double_t ca = TMath::Cos(alpha);
1379 Double_t sa = TMath::Sin(alpha);
1380 Double_t r1 = fC*fX - fE;
46d29e70 1381
15ed8ba1 1382 fX = x1*ca + y1*sa;
1383 fY = -x1*sa + y1*ca;
1384 if ((r1*r1) > 1.0) {
46d29e70 1385 return 0;
1386 }
15ed8ba1 1387 fE = fE*ca + (fC*y1 + TMath::Sqrt(1.0 - r1*r1)) * sa;
46d29e70 1388
15ed8ba1 1389 Double_t r2 = fC*fX - fE;
1390 if (TMath::Abs(r2) >= 0.9) {
1391 Int_t n = GetNumberOfClusters();
1392 if (n > 4) {
1393 AliError(Form("Rotation failed N = %d !\n",n));
1394 }
1395 return 0;
1396 }
46d29e70 1397
15ed8ba1 1398 if ((r2*r2) > 1.0) {
1399 return 0;
1400 }
1401 Double_t y0 = fY + TMath::Sqrt(1.0 - r2*r2) / fC;
1402 if ((fY-y0)*fC >= 0.0) {
1403 Int_t n = GetNumberOfClusters();
1404 if (n > 4) {
1405 AliError(Form("Rotation failed N = %d !\n",n));
1406 }
1407 return 0;
1408 }
46d29e70 1409
15ed8ba1 1410 // f = F - 1
1411 Double_t f00 = ca-1.0;
1412 Double_t f24 = (y1 - r1*x1/TMath::Sqrt(1.0 - r1*r1)) * sa;
1413 Double_t f20 = fC*sa;
1414 Double_t f22 = (ca + sa*r1/TMath::Sqrt(1.0 - r1*r1)) - 1.0;
1415
1416 // b = C*ft
1417 Double_t b00 = fCyy*f00;
1418 Double_t b02 = fCyy*f20 + fCcy*f24 + fCey*f22;
1419 Double_t b10 = fCzy*f00;
1420 Double_t b12 = fCzy*f20 + fCcz*f24 + fCez*f22;
1421 Double_t b20 = fCey*f00;
1422 Double_t b22 = fCey*f20 + fCce*f24 + fCee*f22;
1423 Double_t b30 = fCty*f00;
1424 Double_t b32 = fCty*f20 + fCct*f24 + fCte*f22;
1425 Double_t b40 = fCcy*f00;
1426 Double_t b42 = fCcy*f20 + fCcc*f24 + fCce*f22;
1427
1428 // a = f*b = f*C*ft
1429 Double_t a00 = f00*b00;
1430 Double_t a02 = f00*b02;
1431 Double_t a22 = f20*b02 + f24*b42 + f22*b22;
1432
1433 // F*C*Ft = C + (a + b + bt)
1434 fCyy += a00 + 2.0*b00;
46d29e70 1435 fCzy += b10;
15ed8ba1 1436 fCey += a02 + b20 + b02;
b3a5a838 1437 fCty += b30;
1438 fCcy += b40;
1439 fCez += b12;
1440 fCte += b32;
15ed8ba1 1441 fCee += a22 + 2.0*b22;
b3a5a838 1442 fCce += b42;
46d29e70 1443
b8dc2353 1444 return 1;
46d29e70 1445
53c17fbf 1446}
7ad19338 1447
46d29e70 1448//_____________________________________________________________________________
fd621f36 1449Double_t AliTRDtrack::GetPredictedChi2(const AliTRDcluster *c, Double_t h01) const
46d29e70 1450{
53c17fbf 1451 //
1452 // Returns the track chi2
1453 //
1454
b8dc2353 1455 Bool_t fNoTilt = kTRUE;
15ed8ba1 1456 if (TMath::Abs(h01) > 0.003) {
1457 fNoTilt = kFALSE;
1458 }
1459
1460 Double_t chi2;
1461 Double_t dy;
1462 Double_t r00;
1463 Double_t r01;
1464 Double_t r11;
1465
1466 if (fNoTilt) {
1467
1468 dy = c->GetY() - fY;
1469 r00 = c->GetSigmaY2();
1470 chi2 = (dy*dy) / r00;
b8dc2353 1471
46d29e70 1472 }
b8dc2353 1473 else {
b8dc2353 1474
15ed8ba1 1475 Double_t padlength = TMath::Sqrt(c->GetSigmaZ2() * 12.0);
1476
1477 r00 = c->GetSigmaY2();
1478 r01 = 0.0;
1479 r11 = c->GetSigmaZ2();
1480 r00 += fCyy;
1481 r01 += fCzy;
1482 r11 += fCzz;
1483
1484 Double_t det = r00*r11 - r01*r01;
b8dc2353 1485 if (TMath::Abs(det) < 1.e-10) {
15ed8ba1 1486 Int_t n = GetNumberOfClusters();
1487 if (n > 4) {
1488 AliError(Form("Singular matrix N = %d!\n",n));
1489 }
b8dc2353 1490 return 1e10;
1491 }
15ed8ba1 1492
1493 Double_t tmp = r00;
1494 r00 = r11;
1495 r11 = tmp;
1496 r01 = -r01;
1497 Double_t dy = c->GetY() - fY;
1498 Double_t dz = c->GetZ() - fZ;
4f1c04d3 1499 Double_t tiltdz = dz;
15ed8ba1 1500 if (TMath::Abs(tiltdz) > padlength/2.0) {
1501 tiltdz = TMath::Sign(padlength/2.0,dz);
4f1c04d3 1502 }
15ed8ba1 1503 dy = dy + h01*tiltdz;
1504
1505 chi2 = (dy*r00*dy + 2*r01*dy*dz + dz*r11*dz) / det;
a819a5f7 1506
b8dc2353 1507 }
53c17fbf 1508
b8dc2353 1509 return chi2;
46d29e70 1510
53c17fbf 1511}
7ad19338 1512
46d29e70 1513//_________________________________________________________________________
1514void AliTRDtrack::GetPxPyPz(Double_t& px, Double_t& py, Double_t& pz) const
1515{
15ed8ba1 1516 //
46d29e70 1517 // Returns reconstructed track momentum in the global system.
15ed8ba1 1518 //
46d29e70 1519
15ed8ba1 1520 Double_t pt = TMath::Abs(GetPt());
1521 Double_t r = fC*fX - fE;
5443e65e 1522
1523 Double_t y0;
15ed8ba1 1524 if (r > 1) {
1525 py = pt;
1526 px = 0.0;
1527 }
1528 else if (r < -1) {
1529 py = -pt;
1530 px = 0.0;
1531 }
5443e65e 1532 else {
15ed8ba1 1533 y0 = fY + TMath::Sqrt(1.0 - r*r) / fC;
1534 px = -pt * (fY - y0) * fC; //cos(phi);
1535 py = -pt * (fE - fX*fC); //sin(phi);
5443e65e 1536 }
15ed8ba1 1537
1538 pz = pt*fT;
1539 Double_t tmp = px * TMath::Cos(fAlpha)
1540 - py * TMath::Sin(fAlpha);
1541 py = px * TMath::Sin(fAlpha)
1542 + py * TMath::Cos(fAlpha);
1543 px = tmp;
46d29e70 1544
1545}
1546
5443e65e 1547//_________________________________________________________________________
1548void AliTRDtrack::GetGlobalXYZ(Double_t& x, Double_t& y, Double_t& z) const
46d29e70 1549{
15ed8ba1 1550 //
5443e65e 1551 // Returns reconstructed track coordinates in the global system.
15ed8ba1 1552 //
5443e65e 1553
15ed8ba1 1554 x = fX;
1555 y = fY;
1556 z = fZ
1557;
1558 Double_t tmp = x * TMath::Cos(fAlpha)
1559 - y * TMath::Sin(fAlpha);
1560 y = x * TMath::Sin(fAlpha)
1561 + y * TMath::Cos(fAlpha);
1562 x = tmp;
5443e65e 1563
1564}
3ab6f951 1565
5443e65e 1566//_________________________________________________________________________
53c17fbf 1567void AliTRDtrack::ResetCovariance()
1568{
5443e65e 1569 //
1570 // Resets covariance matrix
1571 //
46d29e70 1572
15ed8ba1 1573 fCyy *= 10.0;
1574 fCzy = 0.0; fCzz *= 10.0;
1575 fCey = 0.0; fCez = 0.0; fCee *= 10.0;
1576 fCty = 0.0; fCtz = 0.0; fCte = 0.0; fCtt *= 10.0;
1577 fCcy = 0.0; fCcz = 0.0; fCce = 0.0; fCct = 0.0; fCcc *= 10.0;
53c17fbf 1578
5443e65e 1579}
b8dc2353 1580
53c17fbf 1581//_____________________________________________________________________________
1582void AliTRDtrack::ResetCovariance(Float_t mult)
1583{
46e2d86c 1584 //
1585 // Resets covariance matrix
1586 //
1587
15ed8ba1 1588 fCyy *= mult;
1589 fCzy *= 0.0; fCzz *= 1.0;
1590 fCey *= 0.0; fCez *= 0.0; fCee *= mult;
1591 fCty *= 0.0; fCtz *= 0.0; fCte *= 0.0; fCtt *= 1.0;
1592 fCcy *= 0.0; fCcz *= 0.0; fCce *= 0.0; fCct *= 0.0; fCcc *= mult;
7ad19338 1593
53c17fbf 1594}
7ad19338 1595
53c17fbf 1596//_____________________________________________________________________________
16d9fbba 1597void AliTRDtrack::MakeBackupTrack()
1598{
1599 //
53c17fbf 1600 // Creates a backup track
16d9fbba 1601 //
53c17fbf 1602
15ed8ba1 1603 if (fBackupTrack) {
1604 delete fBackupTrack;
1605 }
1606
16d9fbba 1607 fBackupTrack = new AliTRDtrack(*this);
4f1c04d3 1608
1609}
1610
53c17fbf 1611//_____________________________________________________________________________
1612Int_t AliTRDtrack::GetProlongation(Double_t xk, Double_t &y, Double_t &z)
1613{
4f1c04d3 1614 //
1615 // Find prolongation at given x
15ed8ba1 1616 // Return 0 if it does not exist
1617 //
4f1c04d3 1618
15ed8ba1 1619 Double_t c1 = fC*fX - fE;
1620 if (TMath::Abs(c1) > 1.0) {
1621 return 0;
1622 }
1623
1624 Double_t r1 = TMath::Sqrt(1.0 - c1*c1);
1625 Double_t c2 = fC*xk - fE;
1626 if (TMath::Abs(c2) > 1.0) {
1627 return 0;
1628 }
1629
1630 Double_t r2 = TMath::Sqrt(1.0 - c2*c2);
1631 y = fY + (xk-fX)*(c1+c2)/(r1+r2);
1632 z = fZ + (xk-fX)*(c1+c2)/(c1*r2 + c2*r1)*fT;
4f1c04d3 1633
1634 return 1;
1635
16d9fbba 1636}
3fad3d32 1637
53c17fbf 1638//_____________________________________________________________________________
15ed8ba1 1639Int_t AliTRDtrack::PropagateToX(Double_t xr, Double_t step)
3fad3d32 1640{
1641 //
15ed8ba1 1642 // Propagate track to a given x position
1643 // Works inside of the 20 degree segmentation
1644 // (local cooordinate frame for TRD , TPC, TOF)
3fad3d32 1645 //
15ed8ba1 1646 // The material budget is taken from the geo manager
3fad3d32 1647 //
15ed8ba1 1648
1649 Double_t xyz0[3];
1650 Double_t xyz1[3];
1651 Double_t y;
1652 Double_t z;
1653
1654 // Critical alpha - cross sector indication
1655 const Double_t kAlphac = TMath::Pi()/9.0;
1656 const Double_t kTalphac = TMath::Tan(kAlphac*0.5);
1657
1658 // Direction +-
1659 Double_t dir = (fX > xr) ? -1.0 : 1.0;
1660
1661 for (Double_t x = fX + dir*step; dir*x < dir*xr; x += dir*step) {
1662
3fad3d32 1663 GetGlobalXYZ(xyz0[0],xyz0[1],xyz0[2]);
1664 GetProlongation(x,y,z);
15ed8ba1 1665 xyz1[0] = x * TMath::Cos(fAlpha) + y * TMath::Sin(fAlpha);
1666 xyz1[1] = x * TMath::Sin(fAlpha) - y * TMath::Cos(fAlpha);
3fad3d32 1667 xyz1[2] = z;
1668 Double_t param[7];
1669 AliKalmanTrack::MeanMaterialBudget(xyz0,xyz1,param);
15ed8ba1 1670
1671 if ((param[0] > 0) &&
1672 (param[1] > 0)) {
1673 PropagateTo(x,param[1],param[0]);
1674 }
1675 if (fY > fX*kTalphac) {
53c17fbf 1676 Rotate(-kAlphac);
3fad3d32 1677 }
15ed8ba1 1678 if (fY < -fX*kTalphac) {
53c17fbf 1679 Rotate(kAlphac);
3fad3d32 1680 }
15ed8ba1 1681
3fad3d32 1682 }
15ed8ba1 1683
3fad3d32 1684 PropagateTo(xr);
53c17fbf 1685
3fad3d32 1686 return 0;
3fad3d32 1687
53c17fbf 1688}
3fad3d32 1689
53c17fbf 1690//_____________________________________________________________________________
15ed8ba1 1691Int_t AliTRDtrack::PropagateToR(Double_t r,Double_t step)
3fad3d32 1692{
1693 //
15ed8ba1 1694 // Propagate a track to a given radial position
1695 // The rotation is always connected to the last track position
3fad3d32 1696 //
15ed8ba1 1697
1698 Double_t xyz0[3];
1699 Double_t xyz1[3];
1700 Double_t y;
1701 Double_t z;
1702
1703 // Direction +-
1704 Double_t radius = TMath::Sqrt(fX*fX + fY*fY);
1705 Double_t dir = (radius > r) ? -1.0 : 1.0;
1706
1707 for (Double_t x = radius + dir*step; dir*x < dir*r; x += dir*step) {
3fad3d32 1708 GetGlobalXYZ(xyz0[0],xyz0[1],xyz0[2]);
1709 Double_t alpha = TMath::ATan2(xyz0[1],xyz0[0]);
1710 Rotate(alpha,kTRUE);
1711 GetGlobalXYZ(xyz0[0],xyz0[1],xyz0[2]);
1712 GetProlongation(x,y,z);
15ed8ba1 1713 xyz1[0] = x * TMath::Cos(alpha) + y * TMath::Sin(alpha);
1714 xyz1[1] = x * TMath::Sin(alpha) - y * TMath::Cos(alpha);
3fad3d32 1715 xyz1[2] = z;
1716 Double_t param[7];
1717 AliKalmanTrack::MeanMaterialBudget(xyz0,xyz1,param);
15ed8ba1 1718 if (param[1] <= 0.0) {
1719 param[1] = 100000000.0;
1720 }
3fad3d32 1721 PropagateTo(x,param[1],param[0]);
1722 }
15ed8ba1 1723
3fad3d32 1724 GetGlobalXYZ(xyz0[0],xyz0[1],xyz0[2]);
1725 Double_t alpha = TMath::ATan2(xyz0[1],xyz0[0]);
1726 Rotate(alpha,kTRUE);
1727 GetGlobalXYZ(xyz0[0],xyz0[1],xyz0[2]);
1728 GetProlongation(r,y,z);
15ed8ba1 1729 xyz1[0] = r * TMath::Cos(alpha) + y * TMath::Sin(alpha);
1730 xyz1[1] = r * TMath::Sin(alpha) - y * TMath::Cos(alpha);
3fad3d32 1731 xyz1[2] = z;
1732 Double_t param[7];
1733 AliKalmanTrack::MeanMaterialBudget(xyz0,xyz1,param);
15ed8ba1 1734
1735 if (param[1] <= 0.0) {
1736 param[1] = 100000000.0;
1737 }
3fad3d32 1738 PropagateTo(r,param[1],param[0]);
53c17fbf 1739
3fad3d32 1740 return 0;
53c17fbf 1741
1742}
1743
1744//_____________________________________________________________________________
4e28c495 1745Int_t AliTRDtrack::GetSector() const
53c17fbf 1746{
1747 //
1748 // Return the current sector
1749 //
1750
1751 return Int_t(TVector2::Phi_0_2pi(fAlpha)
1752 / AliTRDgeometry::GetAlpha())
1753 % AliTRDgeometry::kNsect;
1754
3fad3d32 1755}
1756
53c17fbf 1757//_____________________________________________________________________________
15ed8ba1 1758Double_t AliTRDtrack::Get1Pt() const
53c17fbf 1759{
15ed8ba1 1760 //
5773defd 1761 // Returns the inverse Pt (1/GeV/c)
1762 // (or 1/"most probable pt", if the field is too weak)
15ed8ba1 1763 //
1764
1765 if (TMath::Abs(GetLocalConvConst()) > kVeryBigConvConst) {
1766 return 1.0 / kMostProbableMomentum
1767 / TMath::Sqrt(1.0 + GetTgl()*GetTgl());
1768 }
1769
1770 return (TMath::Sign(1.0e-9,fC) + fC) * GetLocalConvConst();
1771
53c17fbf 1772}
1773
1774//_____________________________________________________________________________
15ed8ba1 1775Double_t AliTRDtrack::GetP() const
53c17fbf 1776{
1777 //
1778 // Returns the total momentum
1779 //
1780
15ed8ba1 1781 return TMath::Abs(GetPt()) * TMath::Sqrt(1.0 + GetTgl()*GetTgl());
53c17fbf 1782
1783}
1784
1785//_____________________________________________________________________________
4e28c495 1786Double_t AliTRDtrack::GetYat(Double_t xk) const
53c17fbf 1787{
1788 //
1789 // This function calculates the Y-coordinate of a track at
1790 // the plane x = xk.
1791 // Needed for matching with the TOF (I.Belikov)
1792 //
1793
1794 Double_t c1 = fC*fX - fE;
1795 Double_t r1 = TMath::Sqrt(1.0 - c1*c1);
1796 Double_t c2 = fC*xk - fE;
1797 Double_t r2 = TMath::Sqrt(1.0- c2*c2);
15ed8ba1 1798
53c17fbf 1799 return fY + (xk-fX)*(c1+c2)/(r1+r2);
1800
1801}
1802
1803//_____________________________________________________________________________
4e28c495 1804void AliTRDtrack::SetSampledEdx(Float_t q, Int_t i)
53c17fbf 1805{
1806 //
1807 // The sampled energy loss
1808 //
1809
1810 Double_t s = GetSnp();
1811 Double_t t = GetTgl();
15ed8ba1 1812
1813 q *= TMath::Sqrt((1.0 - s*s) / (1.0 + t*t));
53c17fbf 1814 fdQdl[i] = q;
1815
1816}
1817
1818 //_____________________________________________________________________________
4e28c495 1819void AliTRDtrack::SetSampledEdx(Float_t q)
53c17fbf 1820{
1821 //
1822 // The sampled energy loss
1823 //
1824
1825 Double_t s = GetSnp();
1826 Double_t t = GetTgl();
15ed8ba1 1827
1828 q*= TMath::Sqrt((1.0 - s*s) / (1.0 + t*t));
53c17fbf 1829 fdQdl[fNdedx] = q;
1830 fNdedx++;
1831
1832}
1833
1834//_____________________________________________________________________________
4e28c495 1835void AliTRDtrack::GetXYZ(Float_t r[3]) const
53c17fbf 1836{
15ed8ba1 1837 //
53c17fbf 1838 // Returns the position of the track in the global coord. system
15ed8ba1 1839 //
53c17fbf 1840
1841 Double_t cs = TMath::Cos(fAlpha);
1842 Double_t sn = TMath::Sin(fAlpha);
1843 r[0] = fX*cs - fY*sn;
1844 r[1] = fX*sn + fY*cs;
1845 r[2] = fZ;
1846
1847}