]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TRD/AliTRDtrack.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Corrected delete statement (Christian)
[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 * *
43bfc8af 7 * Permission to use , copy, modify and distribute this software and its *
46d29e70 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>
fe33d239 19
53c17fbf 20#include <TMath.h>
e1e6896f 21#include <TVector2.h>
46d29e70 22
6c94f330 23#include "AliTracker.h"
53c17fbf 24#include "AliESDtrack.h"
fe33d239 25
46d29e70 26#include "AliTRDgeometry.h"
27#include "AliTRDcluster.h"
28#include "AliTRDtrack.h"
53c17fbf 29#include "AliTRDtracklet.h"
b3a5a838 30
43bfc8af 31// A. Bercuci - used for PID calculations
32#include "AliTRDcalibDB.h"
720a0a16 33#include "Cal/AliTRDCalPID.h"
43bfc8af 34
46d29e70 35ClassImp(AliTRDtrack)
36
53c17fbf 37///////////////////////////////////////////////////////////////////////////////
38// //
39// Represents a reconstructed TRD track //
40// Local TRD Kalman track //
41// //
42///////////////////////////////////////////////////////////////////////////////
7ad19338 43
fe33d239 44//_____________________________________________________________________________
45AliTRDtrack::AliTRDtrack()
46 :AliKalmanTrack()
47 ,fSeedLab(-1)
48 ,fdEdx(0)
49 ,fDE(0)
79143cbf 50 ,fClusterOwner(kFALSE)
fe33d239 51 ,fStopped(kFALSE)
52 ,fLhElectron(0)
53 ,fNWrong(0)
54 ,fNRotate(0)
55 ,fNCross(0)
56 ,fNExpected(0)
57 ,fNLast(0)
58 ,fNExpectedLast(0)
59 ,fNdedx(0)
60 ,fChi2Last(1e10)
61 ,fBackupTrack(0x0)
f146da82 62{
fe33d239 63 //
64 // AliTRDtrack default constructor
65 //
66
67 for (Int_t i = 0; i < kNplane; i++) {
68 for (Int_t j = 0; j < kNslice; j++) {
69 fdEdxPlane[i][j] = 0.0;
6d45eaef 70 }
f146da82 71 fTimBinPlane[i] = -1;
43bfc8af 72 // A.Bercuci additions
73 fMom[i] = -1.;
74 fSnp[i] = 0.;
75 fTgl[i] = 0.;
76 }
fe33d239 77
78 for (UInt_t i = 0; i < kMAXCLUSTERSPERTRACK; i++) {
79 fIndex[i] = 0;
6c94f330 80 fIndexBackup[i] = 0;
fe33d239 81 fdQdl[i] = 0;
43bfc8af 82 //A.Bercuci additions
83 fClusters[i] = 0x0;
84 }
fe33d239 85
86 for (Int_t i = 0; i < 3; i++) {
87 fBudget[i] = 0;
f146da82 88 }
fe33d239 89
f146da82 90}
6d45eaef 91
46d29e70 92//_____________________________________________________________________________
43bfc8af 93AliTRDtrack::AliTRDtrack(AliTRDcluster *c, Int_t index
fe33d239 94 , const Double_t p[5], const Double_t cov[15]
95 , Double_t x, Double_t alpha)
96 :AliKalmanTrack()
97 ,fSeedLab(-1)
98 ,fdEdx(0)
99 ,fDE(0)
79143cbf 100 ,fClusterOwner(kFALSE)
fe33d239 101 ,fStopped(kFALSE)
102 ,fLhElectron(0)
103 ,fNWrong(0)
104 ,fNRotate(0)
105 ,fNCross(0)
106 ,fNExpected(0)
107 ,fNLast(0)
108 ,fNExpectedLast(0)
109 ,fNdedx(0)
110 ,fChi2Last(1e10)
43bfc8af 111 ,fBackupTrack(0x0)
15ed8ba1 112{
fe33d239 113 //
114 // The main AliTRDtrack constructor.
115 //
116
117 Double_t cnv = 1.0/(GetBz() * kB2C);
118
119 Double_t pp[5] = { p[0]
120 , p[1]
121 , x*p[4] - p[2]
122 , p[3]
123 , p[4]*cnv };
6c94f330 124
125 Double_t c22 = x*x*cov[14] - 2*x*cov[12] + cov[5];
126 Double_t c32 = x*cov[13] - cov[8];
fe33d239 127 Double_t c20 = x*cov[10] - cov[3];
128 Double_t c21 = x*cov[11] - cov[4];
129 Double_t c42 = x*cov[14] - cov[12];
6c94f330 130
fe33d239 131 Double_t cc[15] = { cov[0 ]
132 , cov[1 ], cov[2 ]
133 , c20, c21, c22
134 , cov[6 ], cov[7 ], c32, cov[9 ]
135 , cov[10]*cnv, cov[11]*cnv, c42*cnv, cov[13]*cnv, cov[14]*cnv*cnv };
6c94f330 136
137 Set(x,alpha,pp,cc);
5443e65e 138 SetNumberOfClusters(1);
fe33d239 139 fIndex[0] = index;
43bfc8af 140 fClusters[0] = c; // A.Bercuci additions
5443e65e 141
fe33d239 142 for (Int_t i = 0; i < kNplane; i++) {
143 for (Int_t j = 0; j < kNslice; j++) {
144 fdEdxPlane[i][j] = 0.0;
6d45eaef 145 }
146 fTimBinPlane[i] = -1;
43bfc8af 147 // A.Bercuci additions
148 fMom[i] = -1.;
149 fSnp[i] = 0.;
150 fTgl[i] = 0.;
eab5961e 151 }
a2b90f83 152
5443e65e 153 Double_t q = TMath::Abs(c->GetQ());
fe33d239 154 Double_t s = GetSnp();
155 Double_t t = GetTgl();
156 if (s*s < 1) {
157 q *= TMath::Sqrt((1-s*s)/(1+t*t));
158 }
53c17fbf 159
43bfc8af 160 fdQdl[0] = q;
fe33d239 161 for (UInt_t i = 1; i < kMAXCLUSTERSPERTRACK; i++) {
162 fdQdl[i] = 0;
163 fIndex[i] = 0;
164 fIndexBackup[i] = 0;
43bfc8af 165 // A.Bercuci additions
166 fClusters[i] = 0x0;
15ed8ba1 167 }
fe33d239 168
169 for (Int_t i = 0; i < 3;i++) {
170 fBudget[i] = 0;
171 }
172
46d29e70 173}
174
175//_____________________________________________________________________________
43bfc8af 176AliTRDtrack::AliTRDtrack(const AliTRDtrack &t/*, const Bool_t owner*/)
fe33d239 177 :AliKalmanTrack(t)
178 ,fSeedLab(t.GetSeedLabel())
179 ,fdEdx(t.fdEdx)
180 ,fDE(t.fDE)
79143cbf 181 ,fClusterOwner(kTRUE)
fe33d239 182 ,fStopped(t.fStopped)
183 ,fLhElectron(0)
184 ,fNWrong(t.fNWrong)
185 ,fNRotate(t.fNRotate)
186 ,fNCross(t.fNCross)
187 ,fNExpected(t.fNExpected)
188 ,fNLast(t.fNLast)
189 ,fNExpectedLast(t.fNExpectedLast)
190 ,fNdedx(t.fNdedx)
191 ,fChi2Last(t.fChi2Last)
43bfc8af 192 ,fBackupTrack(0x0)
53c17fbf 193{
46d29e70 194 //
195 // Copy constructor.
196 //
fe33d239 197
198 for (Int_t i = 0; i < kNplane; i++) {
199 for (Int_t j = 0; j < kNslice; j++) {
6d45eaef 200 fdEdxPlane[i][j] = t.fdEdxPlane[i][j];
201 }
202 fTimBinPlane[i] = t.fTimBinPlane[i];
203 fTracklets[i] = t.fTracklets[i];
43bfc8af 204 // A.Bercuci additions
205 fMom[i] = t.fMom[i];
206 fSnp[i] = t.fSnp[i];
207 fTgl[i] = t.fTgl[i];
eab5961e 208 }
46d29e70 209
fe33d239 210 Int_t n = t.GetNumberOfClusters();
6c94f330 211 SetNumberOfClusters(n);
fe33d239 212
213 for (Int_t i = 0; i < n; i++) {
214 fIndex[i] = t.fIndex[i];
215 fIndexBackup[i] = t.fIndex[i];
216 fdQdl[i] = t.fdQdl[i];
43bfc8af 217 // A.Bercuci additions
218 if(fClusterOwner && t.fClusters[i]) fClusters[i] = new AliTRDcluster(*(t.fClusters[i]));
219 else fClusters[i] = t.fClusters[i];
220 }
fe33d239 221
222 for (UInt_t i = n; i < kMAXCLUSTERSPERTRACK; i++) {
223 fdQdl[i] = 0;
224 fIndex[i] = 0;
225 fIndexBackup[i] = 0;
43bfc8af 226 // A.Bercuci additions
227 fClusters[i] = 0x0;
03b0452e 228 }
15ed8ba1 229
fe33d239 230 for (Int_t i = 0; i < 3;i++) {
231 fBudget[i] = t.fBudget[i];
15ed8ba1 232 }
fe33d239 233
5443e65e 234}
235
236//_____________________________________________________________________________
fe33d239 237AliTRDtrack::AliTRDtrack(const AliKalmanTrack &t, Double_t /*alpha*/)
238 :AliKalmanTrack(t)
239 ,fSeedLab(-1)
240 ,fdEdx(t.GetPIDsignal())
241 ,fDE(0)
79143cbf 242 ,fClusterOwner(kFALSE)
fe33d239 243 ,fStopped(kFALSE)
244 ,fLhElectron(0.0)
245 ,fNWrong(0)
246 ,fNRotate(0)
247 ,fNCross(0)
248 ,fNExpected(0)
249 ,fNLast(0)
250 ,fNExpectedLast(0)
251 ,fNdedx(0)
252 ,fChi2Last(0.0)
253 ,fBackupTrack(0x0)
53c17fbf 254{
5443e65e 255 //
fe33d239 256 // Constructor from AliTPCtrack or AliITStrack
5443e65e 257 //
258
6c94f330 259 SetLabel(t.GetLabel());
fe33d239 260 SetChi2(0.0);
6c94f330 261 SetMass(t.GetMass());
5443e65e 262 SetNumberOfClusters(0);
263
fe33d239 264 for (Int_t i = 0; i < kNplane; i++) {
265 for (Int_t j = 0; j < kNslice; j++) {
6d45eaef 266 fdEdxPlane[i][j] = 0.0;
267 }
eab5961e 268 fTimBinPlane[i] = -1;
43bfc8af 269 // A.Bercuci additions
270 fMom[i] = -1.;
271 fSnp[i] = 0.;
272 fTgl[i] = 0.;
eab5961e 273 }
5443e65e 274
fe33d239 275 for (UInt_t i = 0; i < kMAXCLUSTERSPERTRACK; i++) {
276 fdQdl[i] = 0;
277 fIndex[i] = 0;
278 fIndexBackup[i] = 0;
43bfc8af 279 // A.Bercuci additions
280 fClusters[i] = 0x0;
79e94bf8 281 }
4f1c04d3 282
fe33d239 283 for (Int_t i = 0; i < 3; i++) {
284 fBudget[i] = 0;
285 }
286
79e94bf8 287}
53c17fbf 288
79e94bf8 289//_____________________________________________________________________________
fe33d239 290AliTRDtrack::AliTRDtrack(const AliESDtrack &t)
291 :AliKalmanTrack()
292 ,fSeedLab(-1)
293 ,fdEdx(t.GetTRDsignal())
294 ,fDE(0)
79143cbf 295 ,fClusterOwner(kFALSE)
fe33d239 296 ,fStopped(kFALSE)
297 ,fLhElectron(0)
298 ,fNWrong(0)
299 ,fNRotate(0)
300 ,fNCross(0)
301 ,fNExpected(0)
302 ,fNLast(0)
303 ,fNExpectedLast(0)
304 ,fNdedx(0)
305 ,fChi2Last(1e10)
306 ,fBackupTrack(0x0)
53c17fbf 307{
79e94bf8 308 //
309 // Constructor from AliESDtrack
310 //
fe33d239 311
79e94bf8 312 SetLabel(t.GetLabel());
fe33d239 313 SetChi2(0.0);
79e94bf8 314 SetMass(t.GetMass());
1e9bb598 315 SetNumberOfClusters(t.GetTRDclusters(fIndex));
15ed8ba1 316
46e2d86c 317 Int_t ncl = t.GetTRDclusters(fIndexBackup);
fe33d239 318 for (UInt_t i = ncl; i < kMAXCLUSTERSPERTRACK; i++) {
319 fIndexBackup[i] = 0;
320 fIndex[i] = 0;
15ed8ba1 321 }
fe33d239 322
323 for (Int_t i = 0; i < kNplane; i++) {
324 for (Int_t j = 0; j < kNslice; j++) {
6d45eaef 325 fdEdxPlane[i][j] = t.GetTRDsignals(i,j);
326 }
eab5961e 327 fTimBinPlane[i] = t.GetTRDTimBin(i);
43bfc8af 328 // A.Bercuci additions
329 fMom[i] = -1.;
330 fSnp[i] = 0.;
331 fTgl[i] = 0.;
eab5961e 332 }
79e94bf8 333
fe33d239 334 const AliExternalTrackParam *par = &t;
335 if (t.GetStatus() & AliESDtrack::kTRDbackup) {
336 par = t.GetOuterParam();
337 if (!par) {
338 AliError("No backup info!");
339 par = &t;
340 }
c5a8e3df 341 }
6c94f330 342 Set(par->GetX(),par->GetAlpha(),par->GetParameter(),par->GetCovariance());
79e94bf8 343
6c94f330 344
fe33d239 345 for (UInt_t i = 0; i < kMAXCLUSTERSPERTRACK; i++) {
43bfc8af 346 fdQdl[i] = 0;
347 // A.Bercuci additions
348 fClusters[i] = 0x0;
349 }
fe33d239 350
351 for (Int_t i = 0; i < 3; i++) {
352 fBudget[i] = 0;
353 }
5443e65e 354
fe33d239 355 if ((t.GetStatus() & AliESDtrack::kTIME) == 0) {
356 return;
357 }
c630aafd 358
c630aafd 359 StartTimeIntegral();
fe33d239 360 Double_t times[10];
361 t.GetIntegratedTimes(times);
362 SetIntegratedTimes(times);
c630aafd 363 SetIntegratedLength(t.GetIntegratedLength());
364
16d9fbba 365}
366
53c17fbf 367//____________________________________________________________________________
368AliTRDtrack::~AliTRDtrack()
3fad3d32 369{
370 //
53c17fbf 371 // Destructor
372 //
3fad3d32 373
fe33d239 374 if (fBackupTrack) {
375 delete fBackupTrack;
376 }
43bfc8af 377 fBackupTrack = 0x0;
378 if (fClusterOwner){
379 UInt_t icluster=0;
380 while(icluster<kMAXCLUSTERSPERTRACK && fClusters[icluster]){
381 delete fClusters[icluster];
382 fClusters[icluster] = 0x0;
383 icluster++;
384 }
385 }
3fad3d32 386
53c17fbf 387}
388
53c17fbf 389//____________________________________________________________________________
390Float_t AliTRDtrack::StatusForTOF()
7ad19338 391{
53c17fbf 392 //
393 // Defines the status of the TOF extrapolation
394 //
03b0452e 395
fe33d239 396 // Definition of res ????
397 Float_t res = (0.2 + 0.8 * (fN / (fNExpected + 5.0)))
398 * (0.4 + 0.6 * fTracklets[5].GetN() / 20.0);
399 res *= (0.25 + 0.8 * 40.0 / (40.0 + fBudget[2]));
03b0452e 400 return res;
401
fe33d239 402 // This part of the function is never reached ????
403 // What defines these parameters ????
404 Int_t status = 0;
405 if (GetNumberOfClusters() < 20) return 0;
406 if ((fN > 110) &&
407 (fChi2/(Float_t(fN)) < 3)) return 3; // Gold
408 if ((fNLast > 30) &&
409 (fChi2Last/(Float_t(fNLast)) < 3)) return 3; // Gold
410 if ((fNLast > 20) &&
411 (fChi2Last/(Float_t(fNLast)) < 2)) return 3; // Gold
412 if ((fNLast/(fNExpectedLast+3.0) > 0.8) &&
413 (fChi2Last/Float_t(fNLast) < 5) &&
414 (fNLast > 20)) return 2; // Silber
415 if ((fNLast > 5) &&
416 (((fNLast+1.0)/(fNExpectedLast+1.0)) > 0.8) &&
417 (fChi2Last/(fNLast-5.0) < 6)) return 1;
53c17fbf 418
4f1c04d3 419 return status;
53c17fbf 420
4f1c04d3 421}
46d29e70 422
423//_____________________________________________________________________________
53c17fbf 424Int_t AliTRDtrack::Compare(const TObject *o) const
425{
426 //
427 // Compares tracks according to their Y2 or curvature
428 //
46d29e70 429
d49e463b 430 AliTRDtrack *t = (AliTRDtrack *) o;
46d29e70 431
d49e463b 432 Double_t co = TMath::Abs(t->GetC());
433 Double_t c = TMath::Abs(GetC());
46d29e70 434
d49e463b 435 if (c > co) {
436 return 1;
437 }
438 else if (c < co) {
439 return -1;
440 }
46d29e70 441 return 0;
53c17fbf 442
46d29e70 443}
444
a819a5f7 445//_____________________________________________________________________________
43bfc8af 446void AliTRDtrack::CookdEdx(Double_t low, Double_t up)
15ed8ba1 447{
448 //
d49e463b 449 // Calculates the truncated dE/dx within the "low" and "up" cuts.
15ed8ba1 450 //
451
d49e463b 452 // Array to sort the dEdx values according to amplitude
43bfc8af 453
454 Float_t sorted[kMAXCLUSTERSPERTRACK];
455 fdEdx = 0.;
456
15ed8ba1 457 // Require at least 10 clusters for a dedx measurement
c6011b06 458 if (fNdedx < 10) return;
15ed8ba1 459
43bfc8af 460
d49e463b 461 // Can fdQdl be negative ????
43bfc8af 462 for (Int_t i = 0; i < fNdedx; i++) sorted[i] = TMath::Abs(fdQdl[i]);
15ed8ba1 463 // Sort the dedx values by amplitude
c6011b06 464 Int_t *index = new Int_t[fNdedx];
465 TMath::Sort(fNdedx, sorted, index, kFALSE);
466
43bfc8af 467 // Sum up the truncated charge between lower and upper bounds
c6011b06 468 Int_t nl = Int_t(low * fNdedx);
469 Int_t nu = Int_t( up * fNdedx);
43bfc8af 470 for (Int_t i = nl; i <= nu; i++) fdEdx += sorted[index[i]];
471 fdEdx /= (nu - nl + 1.0);
472
473 delete[] index;
474}
475
476//_____________________________________________________________________________
477void AliTRDtrack::CookdEdxTimBin()
478{
479 //
480 // Set fdEdxPlane and fTimBinPlane and also get the
481 // Time bin for Max. Cluster
482 //
483 // Authors:
484 // Prashant Shukla (shukla@physi.uni-heidelberg.de)
485 // Alexandru Bercuci (A.Bercuci@gsi.de)
486
487
488 Double_t maxcharge[AliESDtrack::kNPlane]; // max charge in chamber
489 // number of clusters attached to track per chamber and slice
490 Int_t nCluster[AliESDtrack::kNPlane][AliESDtrack::kNSlice];
491 //number of time bins in chamber
492 Int_t ntb = AliTRDcalibDB::Instance()->GetNumberOfTimeBins();
493 Int_t plane; // plane of current cluster
494 Int_t tb; // time bin of current cluster
495 Int_t slice; // curent slice
496 AliTRDcluster *cluster = 0x0; // pointer to current cluster
497
498 // Reset class and local contors/variables
499 for (Int_t iPlane = 0; iPlane < AliESDtrack::kNPlane; iPlane++) {
500 fTimBinPlane[iPlane] = -1;
501 maxcharge[iPlane] = 0.;
502 for (Int_t iSlice = 0; iSlice < AliESDtrack::kNSlice; iSlice++) {
503 fdEdxPlane[iPlane][iSlice] = 0.;
504 nCluster[iPlane][iSlice] = 0;
505 }
15ed8ba1 506 }
c6011b06 507
43bfc8af 508 // start looping over clusters attached to this track
509 for (Int_t iClus = 0; iClus < GetNumberOfClusters(); iClus++) {
510 cluster = fClusters[iClus]; //(AliTRDcluster*)tracker->GetCluster(fIndex[iClus]);
511 if(!cluster) continue;
512
513 // Read info from current cluster
514 plane = AliTRDgeometry::GetPlane(cluster->GetDetector());
515 if (plane < 0 || plane >= AliESDtrack::kNPlane) {
516 AliError(Form("Wrong plane %d", plane));
517 continue;
518 }
519
520 tb = cluster->GetLocalTimeBin();
521 if(tb == 0 || tb >= ntb){
522 AliWarning(Form("time bin 0 or > %d in cluster %d", ntb, iClus));
523 AliInfo(Form("dQ/dl %f", fdQdl[iClus]));
524 continue;
525 }
526
527 slice = tb * AliESDtrack::kNSlice / ntb;
528
529 fdEdxPlane[plane][slice] += fdQdl[iClus];
530 if(fdQdl[iClus] > maxcharge[plane]) {
531 maxcharge[plane] = fdQdl[iClus];
532 fTimBinPlane[plane] = tb;
533 }
534 nCluster[plane][slice]++;
535 } // End of loop over cluster
536
537
538
539 // Normalize fdEdxPlane to number of clusters and set track segments
540 for (Int_t iPlane = 0; iPlane < AliESDtrack::kNPlane; iPlane++) {
541 for (Int_t iSlice = 0; iSlice < AliESDtrack::kNSlice; iSlice++) {
542 if (nCluster[iPlane][iSlice]) fdEdxPlane[iPlane][iSlice] /= nCluster[iPlane][iSlice];
543 }
544 }
545}
546
547
548//_____________________________________________________________________________
549void AliTRDtrack::SetTrackSegmentDirMom(const Int_t plane)
550{
551 if(plane<0 || plane>= kNplane){
552 AliError(Form("Trying to access out of range plane (%d)", plane));
553 return;
554 }
555
556 fSnp[plane] = GetSnp();
557 fTgl[plane] = GetTgl();
558 Double_t p[3]; GetPxPyPz(p);
559 fMom[plane] = TMath::Sqrt(p[0]*p[0] + p[1]*p[1] + p[2]*p[2]);
560}
561
562//_____________________________________________________________________________
563Float_t AliTRDtrack::GetTrackLengthPlane(Int_t plane) const
564{
565 if(plane < 0 || plane >= kNplane) return 0.;
566 return (AliTRDgeometry::AmThick() + AliTRDgeometry::DrThick())/TMath::Sqrt((1.
567- fSnp[plane]*fSnp[plane]) / (1. + fTgl[plane]*fTgl[plane]));
568
569}
570
571
572//_____________________________________________________________________________
573Int_t AliTRDtrack::CookPID(AliESDtrack *esd)
574{
575 //
576 // This function calculates the PID probabilities based on TRD signals
577 //
578 // The method produces probabilities based on the charge
579 // and the position of the maximum time bin in each layer.
580 // The dE/dx information can be used as global charge or 2 to 3
720a0a16 581 // slices. Check AliTRDCalPID and AliTRDCalPIDRefMaker for the actual
43bfc8af 582 // implementation.
583 //
584 // Author
585 // Alex Bercuci (A.Bercuci@gsi.de) 2nd May 2007
586
587
588 AliTRDcalibDB *calibration = AliTRDcalibDB::Instance();
589 if (!calibration) {
590 AliError("No access to calibration data");
591 return -1;
592 }
593
594 // Retrieve the CDB container class with the probability distributions
720a0a16 595 const AliTRDCalPID *pd = calibration->GetPIDLQObject();
43bfc8af 596 if (!pd) {
720a0a16 597 AliError("No access to AliTRDCalPID");
43bfc8af 598 return -1;
599 }
600
601
602 Double_t p0 = 1./AliPID::kSPECIES;
603 if(AliPID::kSPECIES != 5){
604 AliError("Probabilities array defined only for 5 values. Please modify !!");
605 return -1;
606 }
607 Double_t p[] = {p0, p0, p0, p0, p0};
608 Float_t length; // track segment length in chamber
609 Int_t nPlanePID = 0;
43bfc8af 610 // Skip tracks which have no TRD signal at all
611 if (fdEdx == 0.) return -1;
612
613 for (Int_t iPlane = 0; iPlane < AliTRDgeometry::kNplan; iPlane++) {
614
615 length = (AliTRDgeometry::AmThick() + AliTRDgeometry::DrThick())/TMath::Sqrt((1. - fSnp[iPlane]*fSnp[iPlane]) / (1. + fTgl[iPlane]*fTgl[iPlane]));
616
617 // check data
618 if((fdEdxPlane[iPlane][0] + fdEdxPlane[iPlane][1] + fdEdxPlane[iPlane][2]) <= 0.
619 || fTimBinPlane[iPlane] == -1.) continue;
620
43bfc8af 621 // this track segment has fulfilled all requierments
622 nPlanePID++;
623
624 // Get the probabilities for the different particle species
625 for (Int_t iSpecies = 0; iSpecies < AliPID::kSPECIES; iSpecies++) {
626 p[iSpecies] *= pd->GetProbability(iSpecies, fMom[iPlane], fdEdxPlane[iPlane], length);
627 //p[iSpecies] *= pd->GetProbabilityT(iSpecies, fMom[iPlane], fTimBinPlane[iPlane]);
628 }
629 }
630 if(nPlanePID == 0) return 0;
631
632 // normalize probabilities
633 Double_t probTotal = 0.;
634 for (Int_t iSpecies = 0; iSpecies < AliPID::kSPECIES; iSpecies++) probTotal += p[iSpecies];
635
636 if(probTotal <= 0.) {
637 AliWarning("The total probability over all species <= 0. This may be caused by some error in the reference histograms.");
638 return -1;
639 }
640 for(Int_t iSpecies = 0; iSpecies < AliPID::kSPECIES; iSpecies++) p[iSpecies] /= probTotal;
641
642
643 // book PID to the ESD track
644 esd->SetTRDpid(p);
645 esd->SetTRDpidQuality(nPlanePID);
646
647 return 0;
648}
649
3551db50 650
a819a5f7 651
46d29e70 652//_____________________________________________________________________________
fe33d239 653Bool_t AliTRDtrack::PropagateTo(Double_t xk, Double_t x0, Double_t rho)
46d29e70 654{
fe33d239 655 //
46d29e70 656 // Propagates a track of particle with mass=pm to a reference plane
657 // defined by x=xk through media of density=rho and radiationLength=x0
fe33d239 658 //
659
660 if (xk == GetX()) {
661 return kTRUE;
662 }
15ed8ba1 663
fe33d239 664 Double_t oldX = GetX();
665 Double_t oldY = GetY();
666 Double_t oldZ = GetZ();
15ed8ba1 667
fe33d239 668 Double_t bz = GetBz();
46d29e70 669
fe33d239 670 if (!AliExternalTrackParam::PropagateTo(xk,bz)) {
671 return kFALSE;
672 }
9c9d2487 673
fe33d239 674 Double_t x = GetX();
675 Double_t y = GetY();
676 Double_t z = GetZ();
677 Double_t d = TMath::Sqrt((x-oldX)*(x-oldX)
678 + (y-oldY)*(y-oldY)
679 + (z-oldZ)*(z-oldZ));
680 if (oldX < xk) {
681 if (IsStartedTimeIntegral()) {
682 Double_t l2 = d;
683 Double_t crv = GetC();
684 if (TMath::Abs(l2*crv) > 0.0001) {
685 // Make correction for curvature if neccesary
686 l2 = 0.5 * TMath::Sqrt((x-oldX)*(x-oldX) + (y-oldY)*(y-oldY));
687 l2 = 2.0 * TMath::ASin(l2 * crv) / crv;
688 l2 = TMath::Sqrt(l2*l2 + (z-oldZ)*(z-oldZ));
689 }
690 AddTimeStep(l2);
6c94f330 691 }
46d29e70 692 }
15ed8ba1 693
6c94f330 694 Double_t ll = (oldX < xk) ? -d : d;
fe33d239 695 if (!AliExternalTrackParam::CorrectForMaterial(ll*rho/x0,x0,GetMass())) {
696 return kFALSE;
697 }
15ed8ba1 698
fe33d239 699 {
700
701 // Energy losses************************
702 Double_t p2 = (1.0 + GetTgl()*GetTgl()) / (Get1Pt()*Get1Pt());
703 Double_t beta2 = p2 / (p2 + GetMass()*GetMass());
d248a39d 704 if (beta2<1.0e-10 || (5940.0 * beta2/(1.0 - beta2 + 1.0e-10) - beta2) < 0.0) {
fe33d239 705 return kFALSE;
706 }
707
708 Double_t dE = 0.153e-3 / beta2
709 * (log(5940.0 * beta2/(1.0 - beta2 + 1.0e-10)) - beta2)
710 * d * rho;
711 Float_t budget = d * rho;
712 fBudget[0] += budget;
713
714 /*
715 // Suspicious part - think about it ?
716 Double_t kinE = TMath::Sqrt(p2);
717 if (dE > 0.8*kinE) dE = 0.8 * kinE; //
718 if (dE < 0) dE = 0.0; // Not valid region for Bethe bloch
719 */
720
721 fDE += dE;
722
723 /*
724 // Suspicious ! I.B.
725 Double_t sigmade = 0.07 * TMath::Sqrt(TMath::Abs(dE)); // Energy loss fluctuation
726 Double_t sigmac2 = sigmade*sigmade*fC*fC*(p2+GetMass()*GetMass())/(p2*p2);
727 fCcc += sigmac2;
728 fCee += fX*fX * sigmac2;
729 */
6c94f330 730
0d5b5c27 731 }
732
6c94f330 733 return kTRUE;
6d45eaef 734
46d29e70 735}
736
46d29e70 737//_____________________________________________________________________________
fe33d239 738Bool_t AliTRDtrack::Update(const AliTRDcluster *c, Double_t chisq, Int_t index
739 , Double_t h01)
46d29e70 740{
fe33d239 741 //
6c94f330 742 // Assignes found cluster to the track and updates track information
fe33d239 743 //
46d29e70 744
b8dc2353 745 Bool_t fNoTilt = kTRUE;
fe33d239 746 if (TMath::Abs(h01) > 0.003) {
747 fNoTilt = kFALSE;
748 }
749
750 // Add angular effect to the error contribution - MI
6c94f330 751 Float_t tangent2 = GetSnp()*GetSnp();
fe33d239 752 if (tangent2 < 0.90000) {
753 tangent2 = tangent2 / (1.0 - tangent2);
b8dc2353 754 }
fe33d239 755 //Float_t errang = tangent2 * 0.04;
b8dc2353 756
fe33d239 757 Double_t p[2] = {c->GetY(), c->GetZ() };
758 //Double_t cov[3] = {c->GetSigmaY2()+errang, 0.0, c->GetSigmaZ2()*100.0 };
759 Double_t sy2 = c->GetSigmaY2() * 4.0;
760 Double_t sz2 = c->GetSigmaZ2() * 4.0;
761 Double_t cov[3] = {sy2 + h01*h01*sz2, h01*(sy2-sz2), sz2 + h01*h01*sy2 };
46e2d86c 762
fe33d239 763 if (!AliExternalTrackParam::Update(p,cov)) {
764 return kFALSE;
765 }
46d29e70 766
fe33d239 767 Int_t n = GetNumberOfClusters();
768 fIndex[n] = index;
b8dc2353 769 SetNumberOfClusters(n+1);
43bfc8af 770
771
b8dc2353 772 SetChi2(GetChi2()+chisq);
5443e65e 773
6c94f330 774 return kTRUE;
fe33d239 775
776}
53c17fbf 777
46e2d86c 778//_____________________________________________________________________________
43bfc8af 779Int_t AliTRDtrack::UpdateMI(AliTRDcluster *c, Double_t chisq, Int_t index, Double_t h01, Int_t /*plane*/)
fe33d239 780{
781 //
46e2d86c 782 // Assignes found cluster to the track and updates track information
fe33d239 783 //
784
46e2d86c 785 Bool_t fNoTilt = kTRUE;
fe33d239 786 if (TMath::Abs(h01) > 0.003) {
787 fNoTilt = kFALSE;
788 }
789
790 // Add angular effect to the error contribution and make correction - MI
6c94f330 791 Double_t tangent2 = GetSnp()*GetSnp();
792 if (tangent2 < 0.90000){
fe33d239 793 tangent2 = tangent2 / (1.0-tangent2);
15ed8ba1 794 }
6c94f330 795 Double_t tangent = TMath::Sqrt(tangent2);
fe33d239 796 if (GetSnp() < 0) {
797 tangent *= -1;
798 }
799
800 //
801 // Is the following still needed ????
802 //
803
6c94f330 804 // Double_t correction = 0*plane;
805 /*
806 Double_t errang = tangent2*0.04; //
807 Double_t errsys =0.025*0.025*20; //systematic error part
15ed8ba1 808
6c94f330 809 Float_t extend =1;
810 if (c->GetNPads()==4) extend=2;
811 */
812 //if (c->GetNPads()==5) extend=3;
813 //if (c->GetNPads()==6) extend=3;
814 //if (c->GetQ()<15) return 1;
15ed8ba1 815
c5a8e3df 816 /*
6c94f330 817 if (corrector!=0){
3c625a9b 818 //if (0){
819 correction = corrector->GetCorrection(plane,c->GetLocalTimeBin(),tangent);
820 if (TMath::Abs(correction)>0){
821 //if we have info
822 errang = corrector->GetSigma(plane,c->GetLocalTimeBin(),tangent);
823 errang *= errang;
824 errang += tangent2*0.04;
825 }
826 }
c5a8e3df 827 */
6c94f330 828 //
829 //Double_t padlength = TMath::Sqrt(c->GetSigmaZ2()*12.);
830 /*
831 {
832 Double_t dy=c->GetY() - GetY(), dz=c->GetZ() - GetZ();
833 printf("%e %e %e %e\n",dy,dz,padlength/2,h01);
c84a5e9e 834 }
6c94f330 835 */
15ed8ba1 836
fe33d239 837 Double_t p[2] = { c->GetY(), c->GetZ() };
838 //Double_t cov[3]={ (c->GetSigmaY2()+errang+errsys)*extend, 0.0, c->GetSigmaZ2()*10000.0 };
839 Double_t sy2 = c->GetSigmaY2() * 4.0;
840 Double_t sz2 = c->GetSigmaZ2() * 4.0;
841 Double_t cov[3] = { sy2 + h01*h01*sz2, h01*(sy2-sz2), sz2 + h01*h01*sy2 };
46e2d86c 842
fe33d239 843 if (!AliExternalTrackParam::Update(p,cov)) {
844 return kFALSE;
845 }
846
43bfc8af 847 // Register cluster to track
fe33d239 848 Int_t n = GetNumberOfClusters();
849 fIndex[n] = index;
43bfc8af 850 fClusters[n] = c; // A.Bercuci 25.07.07
46e2d86c 851 SetNumberOfClusters(n+1);
fe33d239 852 SetChi2(GetChi2() + chisq);
853
854 return kTRUE;
46e2d86c 855
53c17fbf 856}
7ad19338 857
fe33d239 858// //_____________________________________________________________________________
859// Int_t AliTRDtrack::UpdateMI(const AliTRDtracklet &tracklet)
860// {
861// //
862// // Assignes found tracklet to the track and updates track information
863// //
864// // Can this be removed ????
865// //
866//
867// Double_t r00=(tracklet.GetTrackletSigma2()), r01=0., r11= 10000.;
868// r00+=fCyy; r01+=fCzy; r11+=fCzz;
869// //
870// Double_t det=r00*r11 - r01*r01;
871// Double_t tmp=r00; r00=r11/det; r11=tmp/det; r01=-r01/det;
872// //
873
874// Double_t dy=tracklet.GetY() - fY, dz=tracklet.GetZ() - fZ;
7ad19338 875
15ed8ba1 876
fe33d239 877// Double_t s00 = tracklet.GetTrackletSigma2(); // error pad
878// Double_t s11 = 100000; // error pad-row
879// Float_t h01 = tracklet.GetTilt();
880// //
881// // r00 = fCyy + 2*fCzy*h01 + fCzz*h01*h01+s00;
882// r00 = fCyy + fCzz*h01*h01+s00;
883// // r01 = fCzy + fCzz*h01;
884// r01 = fCzy ;
885// r11 = fCzz + s11;
886// det = r00*r11 - r01*r01;
887// // inverse matrix
888// tmp=r00; r00=r11/det; r11=tmp/det; r01=-r01/det;
889
890// Double_t k00=fCyy*r00+fCzy*r01, k01=fCyy*r01+fCzy*r11;
891// Double_t k10=fCzy*r00+fCzz*r01, k11=fCzy*r01+fCzz*r11;
892// Double_t k20=fCey*r00+fCez*r01, k21=fCey*r01+fCez*r11;
893// Double_t k30=fCty*r00+fCtz*r01, k31=fCty*r01+fCtz*r11;
894// Double_t k40=fCcy*r00+fCcz*r01, k41=fCcy*r01+fCcz*r11;
15ed8ba1 895
fe33d239 896// // K matrix
897// // k00=fCyy*r00+fCzy*(r01+h01*r00),k01=fCyy*r01+fCzy*(r11+h01*r01);
898// // k10=fCzy*r00+fCzz*(r01+h01*r00),k11=fCzy*r01+fCzz*(r11+h01*r01);
899// // k20=fCey*r00+fCez*(r01+h01*r00),k21=fCey*r01+fCez*(r11+h01*r01);
900// // k30=fCty*r00+fCtz*(r01+h01*r00),k31=fCty*r01+fCtz*(r11+h01*r01);
901// // k40=fCcy*r00+fCcz*(r01+h01*r00),k41=fCcy*r01+fCcz*(r11+h01*r01);
902// //
903// //Update measurement
904// Double_t cur=fC + k40*dy + k41*dz, eta=fE + k20*dy + k21*dz;
905// // cur=fC + k40*dy + k41*dz; eta=fE + k20*dy + k21*dz;
906// if (TMath::Abs(cur*fX-eta) >= 0.90000) {
907// //Int_t n=GetNumberOfClusters();
908// // if (n>4) cerr<<n<<" AliTRDtrack warning: Filtering failed !\n";
909// return 0;
910// }
911// // k01+=h01*k00;
912// // k11+=h01*k10;
913// // k21+=h01*k20;
914// // k31+=h01*k30;
915// // k41+=h01*k40;
916
917
918// fY += k00*dy + k01*dz;
919// fZ += k10*dy + k11*dz;
920// fE = eta;
921// fT += k30*dy + k31*dz;
922// fC = cur;
7ad19338 923
6c94f330 924
fe33d239 925// //Update covariance
926// //
927// //
928// Double_t oldyy = fCyy, oldzz = fCzz; //, oldee=fCee, oldcc =fCcc;
929// Double_t oldzy = fCzy, oldey = fCey, oldty=fCty, oldcy =fCcy;
930// Double_t oldez = fCez, oldtz = fCtz, oldcz=fCcz;
931// //Double_t oldte = fCte, oldce = fCce;
932// //Double_t oldct = fCct;
933
934// fCyy-=k00*oldyy+k01*oldzy;
935// fCzy-=k10*oldyy+k11*oldzy;
936// fCey-=k20*oldyy+k21*oldzy;
937// fCty-=k30*oldyy+k31*oldzy;
938// fCcy-=k40*oldyy+k41*oldzy;
939// //
940// fCzz-=k10*oldzy+k11*oldzz;
941// fCez-=k20*oldzy+k21*oldzz;
942// fCtz-=k30*oldzy+k31*oldzz;
943// fCcz-=k40*oldzy+k41*oldzz;
944// //
945// fCee-=k20*oldey+k21*oldez;
946// fCte-=k30*oldey+k31*oldez;
947// fCce-=k40*oldey+k41*oldez;
948// //
949// fCtt-=k30*oldty+k31*oldtz;
950// fCct-=k40*oldty+k41*oldtz;
951// //
952// fCcc-=k40*oldcy+k41*oldcz;
953// //
15ed8ba1 954
fe33d239 955// //Int_t n=GetNumberOfClusters();
956// //fIndex[n]=index;
957// //SetNumberOfClusters(n+1);
6c94f330 958
fe33d239 959// //SetChi2(GetChi2()+chisq);
960// // cerr<<"in update: fIndex["<<fN<<"] = "<<index<<endl;
7ad19338 961
fe33d239 962// return 1;
7ad19338 963
fe33d239 964// }
c84a5e9e 965
46d29e70 966//_____________________________________________________________________________
6c94f330 967Bool_t AliTRDtrack::Rotate(Double_t alpha, Bool_t absolute)
46d29e70 968{
fe33d239 969 //
6c94f330 970 // Rotates track parameters in R*phi plane
971 // if absolute rotation alpha is in global system
972 // otherwise alpha rotation is relative to the current rotation angle
fe33d239 973 //
974
3fad3d32 975 if (absolute) {
6c94f330 976 alpha -= GetAlpha();
3fad3d32 977 }
978 else{
979 fNRotate++;
980 }
46d29e70 981
6c94f330 982 return AliExternalTrackParam::Rotate(GetAlpha()+alpha);
fe33d239 983
53c17fbf 984}
7ad19338 985
46d29e70 986//_____________________________________________________________________________
fd621f36 987Double_t AliTRDtrack::GetPredictedChi2(const AliTRDcluster *c, Double_t h01) const
46d29e70 988{
53c17fbf 989 //
990 // Returns the track chi2
991 //
992
fe33d239 993 Double_t p[2] = { c->GetY(), c->GetZ() };
994 Double_t sy2 = c->GetSigmaY2() * 4.0;
995 Double_t sz2 = c->GetSigmaZ2() * 4.0;
996 Double_t cov[3] = { sy2 + h01*h01*sz2, h01*(sy2-sz2), sz2 + h01*h01*sy2 };
6c94f330 997
998 return AliExternalTrackParam::GetPredictedChi2(p,cov);
15ed8ba1 999
fe33d239 1000 //
1001 // Can the following be removed ????
1002 //
6c94f330 1003 /*
1004 Bool_t fNoTilt = kTRUE;
1005 if(TMath::Abs(h01) > 0.003) fNoTilt = kFALSE;
15ed8ba1 1006
6c94f330 1007 return (c->GetY() - GetY())*(c->GetY() - GetY())/c->GetSigmaY2();
1008 */
15ed8ba1 1009
6c94f330 1010 /*
1011 Double_t chi2, dy, r00, r01, r11;
b8dc2353 1012
6c94f330 1013 if(fNoTilt) {
1014 dy=c->GetY() - fY;
1015 r00=c->GetSigmaY2();
1016 chi2 = (dy*dy)/r00;
46d29e70 1017 }
b8dc2353 1018 else {
6c94f330 1019 Double_t padlength = TMath::Sqrt(c->GetSigmaZ2()*12);
1020 //
1021 r00=c->GetSigmaY2(); r01=0.; r11=c->GetSigmaZ2();
1022 r00+=fCyy; r01+=fCzy; r11+=fCzz;
b8dc2353 1023
6c94f330 1024 Double_t det=r00*r11 - r01*r01;
b8dc2353 1025 if (TMath::Abs(det) < 1.e-10) {
6c94f330 1026 Int_t n=GetNumberOfClusters();
1027 if (n>4) cerr<<n<<" AliTRDtrack warning: Singular matrix !\n";
b8dc2353 1028 return 1e10;
1029 }
6c94f330 1030 Double_t tmp=r00; r00=r11; r11=tmp; r01=-r01;
1031 Double_t dy=c->GetY() - fY, dz=c->GetZ() - fZ;
4f1c04d3 1032 Double_t tiltdz = dz;
6c94f330 1033 if (TMath::Abs(tiltdz)>padlength/2.) {
1034 tiltdz = TMath::Sign(padlength/2,dz);
4f1c04d3 1035 }
6c94f330 1036 // dy=dy+h01*dz;
1037 dy=dy+h01*tiltdz;
a819a5f7 1038
6c94f330 1039 chi2 = (dy*r00*dy + 2*r01*dy*dz + dz*r11*dz)/det;
b8dc2353 1040 }
53c17fbf 1041
b8dc2353 1042 return chi2;
6c94f330 1043 */
fe33d239 1044
1045}
7ad19338 1046
53c17fbf 1047//_____________________________________________________________________________
16d9fbba 1048void AliTRDtrack::MakeBackupTrack()
1049{
1050 //
53c17fbf 1051 // Creates a backup track
16d9fbba 1052 //
53c17fbf 1053
fe33d239 1054 if (fBackupTrack) {
1055 delete fBackupTrack;
1056 }
16d9fbba 1057 fBackupTrack = new AliTRDtrack(*this);
4f1c04d3 1058
1059}
1060
53c17fbf 1061//_____________________________________________________________________________
1062Int_t AliTRDtrack::GetProlongation(Double_t xk, Double_t &y, Double_t &z)
1063{
4f1c04d3 1064 //
fe33d239 1065 // Find a prolongation at given x
1066 // Return 0 if it does not exist
1067 //
1068
1069 Double_t bz = GetBz();
15ed8ba1 1070
fe33d239 1071 if (!AliExternalTrackParam::GetYAt(xk,bz,y)) {
1072 return 0;
1073 }
1074 if (!AliExternalTrackParam::GetZAt(xk,bz,z)) {
1075 return 0;
1076 }
4f1c04d3 1077
6c94f330 1078 return 1;
fe33d239 1079
16d9fbba 1080}
3fad3d32 1081
53c17fbf 1082//_____________________________________________________________________________
fe33d239 1083Int_t AliTRDtrack::PropagateToX(Double_t xr, Double_t step)
3fad3d32 1084{
1085 //
6c94f330 1086 // Propagate track to given x position
fe33d239 1087 // Works inside of the 20 degree segmentation (local cooordinate frame for TRD , TPC, TOF)
3fad3d32 1088 //
fe33d239 1089 // Material budget from geo manager
6c94f330 1090 //
fe33d239 1091
1092 Double_t xyz0[3];
1093 Double_t xyz1[3];
1094 Double_t y;
1095 Double_t z;
1096
1097 const Double_t kAlphac = TMath::Pi()/9.0;
15ed8ba1 1098 const Double_t kTalphac = TMath::Tan(kAlphac*0.5);
fe33d239 1099
1100 // Critical alpha - cross sector indication
1101 Double_t dir = (GetX()>xr) ? -1.0 : 1.0;
1102
1103 // Direction +-
1104 for (Double_t x = GetX()+dir*step; dir*x < dir*xr; x += dir*step) {
1105
6c94f330 1106 GetXYZ(xyz0);
3fad3d32 1107 GetProlongation(x,y,z);
fe33d239 1108 xyz1[0] = x * TMath::Cos(GetAlpha()) + y * TMath::Sin(GetAlpha());
1109 xyz1[1] = x * TMath::Sin(GetAlpha()) - y * TMath::Cos(GetAlpha());
3fad3d32 1110 xyz1[2] = z;
1111 Double_t param[7];
1112 AliKalmanTrack::MeanMaterialBudget(xyz0,xyz1,param);
fe33d239 1113
1114 if ((param[0] > 0) &&
1115 (param[1] > 0)) {
1116 PropagateTo(x,param[1],param[0]);
1117 }
1118
1119 if (GetY() > GetX()*kTalphac) {
53c17fbf 1120 Rotate(-kAlphac);
3fad3d32 1121 }
fe33d239 1122 if (GetY() < -GetX()*kTalphac) {
1123 Rotate( kAlphac);
3fad3d32 1124 }
fe33d239 1125
3fad3d32 1126 }
fe33d239 1127
3fad3d32 1128 PropagateTo(xr);
53c17fbf 1129
3fad3d32 1130 return 0;
3fad3d32 1131
53c17fbf 1132}
3fad3d32 1133
53c17fbf 1134//_____________________________________________________________________________
6c94f330 1135Int_t AliTRDtrack::PropagateToR(Double_t r,Double_t step)
3fad3d32 1136{
1137 //
fe33d239 1138 // Propagate track to the radial position
1139 // Rotation always connected to the last track position
3fad3d32 1140 //
fe33d239 1141
1142 Double_t xyz0[3];
1143 Double_t xyz1[3];
1144 Double_t y;
1145 Double_t z;
1146
6c94f330 1147 Double_t radius = TMath::Sqrt(GetX()*GetX() + GetY()*GetY());
fe33d239 1148 // Direction +-
1149 Double_t dir = (radius>r) ? -1.0 : 1.0;
1150
1151 for (Double_t x = radius+dir*step; dir*x < dir*r; x += dir*step) {
1152
6c94f330 1153 GetXYZ(xyz0);
3fad3d32 1154 Double_t alpha = TMath::ATan2(xyz0[1],xyz0[0]);
1155 Rotate(alpha,kTRUE);
6c94f330 1156 GetXYZ(xyz0);
3fad3d32 1157 GetProlongation(x,y,z);
fe33d239 1158 xyz1[0] = x * TMath::Cos(alpha) + y * TMath::Sin(alpha);
1159 xyz1[1] = x * TMath::Sin(alpha) - y * TMath::Cos(alpha);
3fad3d32 1160 xyz1[2] = z;
1161 Double_t param[7];
1162 AliKalmanTrack::MeanMaterialBudget(xyz0,xyz1,param);
fe33d239 1163 if (param[1] <= 0) {
1164 param[1] = 100000000;
1165 }
3fad3d32 1166 PropagateTo(x,param[1],param[0]);
fe33d239 1167
3fad3d32 1168 }
fe33d239 1169
6c94f330 1170 GetXYZ(xyz0);
3fad3d32 1171 Double_t alpha = TMath::ATan2(xyz0[1],xyz0[0]);
1172 Rotate(alpha,kTRUE);
6c94f330 1173 GetXYZ(xyz0);
3fad3d32 1174 GetProlongation(r,y,z);
fe33d239 1175 xyz1[0] = r * TMath::Cos(alpha) + y * TMath::Sin(alpha);
1176 xyz1[1] = r * TMath::Sin(alpha) - y * TMath::Cos(alpha);
3fad3d32 1177 xyz1[2] = z;
1178 Double_t param[7];
1179 AliKalmanTrack::MeanMaterialBudget(xyz0,xyz1,param);
fe33d239 1180
1181 if (param[1] <= 0) {
1182 param[1] = 100000000;
1183 }
3fad3d32 1184 PropagateTo(r,param[1],param[0]);
53c17fbf 1185
3fad3d32 1186 return 0;
53c17fbf 1187
1188}
1189
1190//_____________________________________________________________________________
4e28c495 1191Int_t AliTRDtrack::GetSector() const
53c17fbf 1192{
1193 //
1194 // Return the current sector
1195 //
1196
fe33d239 1197 return Int_t(TVector2::Phi_0_2pi(GetAlpha()) / AliTRDgeometry::GetAlpha())
1198 % AliTRDgeometry::kNsect;
53c17fbf 1199
3fad3d32 1200}
1201
53c17fbf 1202//_____________________________________________________________________________
4e28c495 1203void AliTRDtrack::SetSampledEdx(Float_t q, Int_t i)
53c17fbf 1204{
1205 //
1206 // The sampled energy loss
1207 //
1208
1209 Double_t s = GetSnp();
1210 Double_t t = GetTgl();
fe33d239 1211 q *= TMath::Sqrt((1.0 - s*s) / (1.0 + t*t));
53c17fbf 1212 fdQdl[i] = q;
1213
1214}
1215
fe33d239 1216//_____________________________________________________________________________
4e28c495 1217void AliTRDtrack::SetSampledEdx(Float_t q)
53c17fbf 1218{
1219 //
1220 // The sampled energy loss
1221 //
1222
1223 Double_t s = GetSnp();
1224 Double_t t = GetTgl();
fe33d239 1225 q *= TMath::Sqrt((1.0 - s*s) / (1.0 + t*t));
53c17fbf 1226 fdQdl[fNdedx] = q;
1227 fNdedx++;
1228
1229}
1230
fe33d239 1231//_____________________________________________________________________________
1232Double_t AliTRDtrack::GetBz() const
1233{
15ed8ba1 1234 //
fe33d239 1235 // Returns Bz component of the magnetic field (kG)
15ed8ba1 1236 //
53c17fbf 1237
fe33d239 1238 if (AliTracker::UniformField()) {
1239 return AliTracker::GetBz();
1240 }
1241 Double_t r[3];
1242 GetXYZ(r);
1243
1244 return AliTracker::GetBz(r);
53c17fbf 1245
fe33d239 1246}