updates in the PWGPP/TRD
[u/mrichter/AliRoot.git] / PWGPP / TRD / info / AliTRDtrackInfo.cxx
CommitLineData
1ee39b3a 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/* $Id: AliTRDtrackInfo.cxx 27496 2008-07-22 08:35:45Z cblume $ */
17
18////////////////////////////////////////////////////////////////////////////
19// //
20// Reconstruction QA //
21// //
22// Authors: //
23// Alex Bercuci <A.Bercuci@gsi.de> //
24// Markus Fasel <M.Fasel@gsi.de> //
25// //
26////////////////////////////////////////////////////////////////////////////
27
28#include "TDatabasePDG.h"
b9ddd472 29#include "TPDGCode.h"
4226db3e 30#include "TVectorT.h"
1ee39b3a 31
32#include "AliTrackReference.h"
eb05d549 33#include "AliTrackPointArray.h"
1ee39b3a 34#include "AliExternalTrackParam.h"
35
36#include "AliTRDseedV1.h"
37#include "AliTRDtrackV1.h"
38#include "AliTRDgeometry.h"
39#include "AliTRDtrackerV1.h"
40
41#include "AliTRDtrackInfo.h"
42
43ClassImp(AliTRDtrackInfo)
44ClassImp(AliTRDtrackInfo::AliMCinfo)
45ClassImp(AliTRDtrackInfo::AliESDinfo)
5066aa9a 46Double_t AliTRDtrackInfo::AliMCinfo::fgKalmanStep = 2.;
3ceb45ae 47Bool_t AliTRDtrackInfo::AliMCinfo::fgKalmanUpdate = kTRUE;
1ee39b3a 48
49//___________________________________________________
50AliTRDtrackInfo::AliTRDtrackInfo():
51 TObject()
52 ,fNClusters(0)
f232df0d 53 ,fTRDtrack(NULL)
54 ,fMC(NULL)
1ee39b3a 55 ,fESD()
56{
57 //
58 // Default constructor
59 //
60}
61
62
63//___________________________________________________
64AliTRDtrackInfo::AliTRDtrackInfo(const AliTRDtrackInfo &trdInfo):
65 TObject((const TObject&)trdInfo)
66 ,fNClusters(trdInfo.fNClusters)
f232df0d 67 ,fTRDtrack(NULL)
a340fe39 68 ,fMC(NULL)
1ee39b3a 69 ,fESD(trdInfo.fESD)
70{
71 //
72 // copy Entries
73 //
74
75 if(trdInfo.fMC) fMC = new AliMCinfo(*trdInfo.fMC);
eb05d549 76 SetTrack(trdInfo.fTRDtrack);
1ee39b3a 77}
78
79//___________________________________________________
80AliTRDtrackInfo::AliMCinfo::AliMCinfo()
81 :fLabel(0)
eb05d549 82 ,fTRDlabel(0)
1ee39b3a 83 ,fPDG(0)
84 ,fNTrackRefs(0)
85{
86 // Set 0-Pointers
87 memset(fTrackRefs, 0, sizeof(AliTrackReference *) * 12);
88}
89
90//___________________________________________________
91AliTRDtrackInfo::AliMCinfo::AliMCinfo(const AliMCinfo &mc)
92 :fLabel(mc.fLabel)
eb05d549 93 ,fTRDlabel(mc.fTRDlabel)
1ee39b3a 94 ,fPDG(mc.fPDG)
95 ,fNTrackRefs(mc.fNTrackRefs)
96{
97 //
98 // Constructor
99 //
100
101 memset(fTrackRefs, 0, sizeof(AliTrackReference *) * 12);
102 for(Int_t ien = 0; ien < 12; ien++){
103 if(mc.fTrackRefs[ien])
104 fTrackRefs[ien] = new AliTrackReference(*(mc.fTrackRefs[ien]));
105 }
106}
107
b9ddd472 108//________________________________________________________
109Int_t AliTRDtrackInfo::AliMCinfo::GetPID() const
110{
61f6b45e 111// Translate pdg code to PID index
112
b9ddd472 113 switch(fPDG){
114 case kElectron:
115 case kPositron: return AliPID::kElectron;
116 case kMuonPlus:
117 case kMuonMinus: return AliPID::kMuon;
118 case kPiPlus:
119 case kPiMinus: return AliPID::kPion;
120 case kKPlus:
121 case kKMinus: return AliPID::kKaon;
122 case kProton:
123 case kProtonBar: return AliPID::kProton;
124 }
125 return -1;
126}
127
1ee39b3a 128//___________________________________________________
129AliTRDtrackInfo::AliESDinfo::AliESDinfo()
b9ddd472 130 :fHasV0(0)
131 ,fId(-1)
1ee39b3a 132 ,fStatus(0)
133 ,fKinkIndex(0)
134 ,fTPCncls(0)
2b33fd9c 135 ,fTPCdedx(0.)
f036f7ec 136 ,fTOFbeta(0.)
3ceb45ae 137 ,fTOFbc(0)
1ee39b3a 138 ,fTRDpidQuality(0)
139 ,fTRDnSlices(0)
eb05d549 140 ,fPt(0.)
141 ,fPhi(-999.)
142 ,fEta(-999.)
f232df0d 143 ,fTRDslices(NULL)
144 ,fOP(NULL)
35983729 145 ,fTPCout(NULL)
b9058d72 146 ,fITSout(NULL)
eb05d549 147 ,fTPArray(NULL)
1ee39b3a 148{
149 //
150 // Constructor
151 //
152
153 memset(fTRDr, 0, AliPID::kSPECIES*sizeof(Double32_t));
d80a6a00 154 memset(fTRDv0pid, 0, AliPID::kSPECIES*sizeof(Int_t));
1ee39b3a 155}
156
157//___________________________________________________
158AliTRDtrackInfo::AliESDinfo::AliESDinfo(const AliESDinfo &esd)
b9ddd472 159 :fHasV0(esd.fHasV0)
160 ,fId(esd.fId)
1ee39b3a 161 ,fStatus(esd.fStatus)
162 ,fKinkIndex(esd.fKinkIndex)
163 ,fTPCncls(esd.fTPCncls)
2b33fd9c 164 ,fTPCdedx(esd.fTPCdedx)
f036f7ec 165 ,fTOFbeta(esd.fTOFbeta)
3ceb45ae 166 ,fTOFbc(esd.fTOFbc)
1ee39b3a 167 ,fTRDpidQuality(esd.fTRDpidQuality)
168 ,fTRDnSlices(esd.fTRDnSlices)
eb05d549 169 ,fPt(esd.fPt)
170 ,fPhi(esd.fPhi)
171 ,fEta(esd.fEta)
f232df0d 172 ,fTRDslices(NULL)
173 ,fOP(NULL)
35983729 174 ,fTPCout(NULL)
b9058d72 175 ,fITSout(NULL)
eb05d549 176 ,fTPArray(NULL)
1ee39b3a 177{
178 //
179 // Constructor
180 //
181
182 memcpy(fTRDr, esd.fTRDr, AliPID::kSPECIES*sizeof(Double32_t));
d80a6a00 183 memcpy(fTRDv0pid, esd.fTRDv0pid, AliPID::kSPECIES*sizeof(Int_t));
1ee39b3a 184
185 if(fTRDnSlices){
186 fTRDslices = new Double32_t[fTRDnSlices];
187 memcpy(fTRDslices, esd.fTRDslices, fTRDnSlices*sizeof(Double32_t));
188 }
eb05d549 189 SetOuterParam(esd.fOP);
190 SetTPCoutParam(esd.fTPCout);
191 SetITSoutParam(esd.fITSout);
192 SetTrackPointArray(esd.fTPArray);
1ee39b3a 193}
194
195
196//___________________________________________________
197AliTRDtrackInfo::~AliTRDtrackInfo()
198{
199 //
200 // Destructor
201 //
202
203 if(fMC) delete fMC;
204 if(fTRDtrack) delete fTRDtrack;
205}
206
207//___________________________________________________
208AliTRDtrackInfo::AliMCinfo::~AliMCinfo()
209{
210 //
211 // Destructor
212 //
213
214 fNTrackRefs = 0;
215 for(Int_t ien = 0; ien < 12; ien++){
216 if(fTrackRefs[ien]) delete fTrackRefs[ien];
f232df0d 217 fTrackRefs[ien] = NULL;
1ee39b3a 218 }
219}
220
221//___________________________________________________
222AliTRDtrackInfo::AliESDinfo::~AliESDinfo()
223{
224 //
225 // Destructor
226 //
227
228 if(fTRDnSlices){
229 delete [] fTRDslices;
f232df0d 230 fTRDslices = NULL;
1ee39b3a 231 fTRDnSlices = 0;
232 }
f232df0d 233 if(fOP) delete fOP; fOP = NULL;
35983729 234 if(fTPCout) delete fTPCout; fTPCout = NULL;
b9058d72 235 if(fITSout) delete fITSout; fITSout = NULL;
eb05d549 236 if(fTPArray) delete fTPArray;
1ee39b3a 237}
238
8778aca3 239//___________________________________________________
240void AliTRDtrackInfo::AliESDinfo::Delete(const Option_t *){
241 //
242 // Delete Pointer members
243 //
244 if(fTRDnSlices){
245 delete [] fTRDslices;
246 fTRDslices = NULL;
247 fTRDnSlices = 0;
248 }
249 if(fOP) delete fOP; fOP = NULL;
35983729 250 if(fTPCout) delete fTPCout; fTPCout = NULL;
b9058d72 251 if(fITSout) delete fITSout; fITSout = NULL;
eb05d549 252 if(fTPArray) delete fTPArray; fTPArray = NULL;
8778aca3 253}
254
1ee39b3a 255
256//___________________________________________________
257AliTRDtrackInfo& AliTRDtrackInfo::operator=(const AliTRDtrackInfo &trdInfo)
258{
259 //
260 // = Operator
261 //
01ccc21a 262 if(this == &trdInfo) return *this;
1ee39b3a 263
264 fNClusters = trdInfo.fNClusters;
1881bb2b 265 fESD = trdInfo.fESD;
1ee39b3a 266
267 if(trdInfo.fMC){
2c33fb46 268 if(!fMC) fMC = new AliMCinfo(*trdInfo.fMC);
269 else{
270 fMC->~AliMCinfo();
1ee39b3a 271 new(fMC) AliMCinfo(*trdInfo.fMC);
2c33fb46 272 }
1881bb2b 273 } else {if(fMC) delete fMC; fMC = NULL;}
1ee39b3a 274
eb05d549 275 SetTrack(trdInfo.fTRDtrack);
1ee39b3a 276
277 return *this;
278}
279
280//___________________________________________________
281AliTRDtrackInfo::AliMCinfo& AliTRDtrackInfo::AliMCinfo::operator=(const AliMCinfo &mc)
282{
283 //
284 // Assignment operator
285 //
286
01ccc21a 287 if(this == &mc) return *this;
1ee39b3a 288 fLabel = mc.fLabel;
eb05d549 289 fTRDlabel = mc.fTRDlabel;
1ee39b3a 290 fPDG = mc.fPDG;
291 fNTrackRefs = mc.fNTrackRefs;
292
293 AliTrackReference **itr = &fTrackRefs[0];
294 AliTrackReference* const *jtr = &mc.fTrackRefs[0];
eb05d549 295 for(Int_t ien = 0; ien < 2*AliTRDgeometry::kNlayer; ien++, itr++, jtr++){
1ee39b3a 296 if((*jtr)){
297 if(!(*itr)) (*itr) = new AliTrackReference(*(*jtr));
2c33fb46 298 else{
299 (*itr)->~AliTrackReference();
300 new(&(*itr)) AliTrackReference(*(*jtr));
301 }
1881bb2b 302 } else {if((*itr)) delete (*itr); (*itr) = NULL;}
1ee39b3a 303 }
304 return *this;
305}
306
307//___________________________________________________
308AliTRDtrackInfo::AliESDinfo& AliTRDtrackInfo::AliESDinfo::operator=(const AliESDinfo &esd)
309{
310 //
311 // Assignment operator
312 //
313
01ccc21a 314 if(this == &esd) return *this;
315 fHasV0 = esd.fHasV0;
1ee39b3a 316 fId = esd.fId;
01ccc21a 317 fStatus = esd.fStatus;
318 fKinkIndex = esd.fKinkIndex;
319 fTPCncls = esd.fTPCncls;
2b33fd9c 320 fTPCdedx = esd.fTPCdedx;
f036f7ec 321 fTOFbeta = esd.fTOFbeta;
01ccc21a 322 fTOFbc = esd.fTOFbc;
323 fTRDpidQuality= esd.fTRDpidQuality;
1ee39b3a 324 fTRDnSlices = esd.fTRDnSlices;
eb05d549 325 fPt = esd.fPt;
326 fPhi = esd.fPhi;
327 fEta = esd.fEta;
178d284a 328
1ee39b3a 329 memcpy(fTRDr, esd.fTRDr, AliPID::kSPECIES*sizeof(Double32_t));
01ccc21a 330 memcpy(fTRDv0pid, esd.fTRDv0pid, AliPID::kSPECIES*sizeof(Int_t));
1ee39b3a 331
332 if(fTRDnSlices){
56838f9d 333 if(!fTRDslices) fTRDslices = new Double32_t[fTRDnSlices];
1ee39b3a 334 memcpy(fTRDslices, esd.fTRDslices, fTRDnSlices*sizeof(Double32_t));
335 }
eb05d549 336 SetOuterParam(esd.fOP);
337 SetTPCoutParam(esd.fTPCout);
338 SetITSoutParam(esd.fITSout);
339 SetTrackPointArray(esd.fTPArray);
1ee39b3a 340
341 return *this;
342}
343
344//___________________________________________________
345void AliTRDtrackInfo::Delete(const Option_t *)
346{
347 //
348 // Delete
349 //
350
2c33fb46 351 AliDebug(2, Form("track[%p] mc[%p]", (void*)fTRDtrack, (void*)fMC));
1ee39b3a 352 fNClusters = 0;
f232df0d 353 if(fMC) delete fMC; fMC = NULL;
8778aca3 354 fESD.Delete(NULL);
f232df0d 355 if(fTRDtrack) delete fTRDtrack; fTRDtrack = NULL;
1ee39b3a 356}
357
358//___________________________________________________
359void AliTRDtrackInfo::SetTrack(const AliTRDtrackV1 *track)
360{
361 //
362 // Set the TRD track
363 //
364
eb05d549 365 if(track){
366 if(!fTRDtrack) fTRDtrack = new AliTRDtrackV1(*track);
367 else{
368 fTRDtrack->~AliTRDtrackV1();
369 new(fTRDtrack) AliTRDtrackV1(*track);
370 }
371 if(track->IsOwner()) fTRDtrack->SetOwner();
372 } else {
373 if(fTRDtrack) delete fTRDtrack; fTRDtrack = NULL;
374 }
375}
376
377//___________________________________________________
378void AliTRDtrackInfo::AliESDinfo::SetTrackPointArray(const AliTrackPointArray *tps)
379{
380 //
381 // Set the track point array for alignment task
382 //
383
384
385 if(tps){
386 if(!fTPArray) fTPArray = new AliTrackPointArray(*tps);
387 else{
388 fTPArray->~AliTrackPointArray();
389 new(fTPArray) AliTrackPointArray(*tps);
390 }
391 } else {
392 if(fTPArray) delete fTPArray; fTPArray = NULL;
2c33fb46 393 }
1ee39b3a 394}
395
396//___________________________________________________
397void AliTRDtrackInfo::AddTrackRef(const AliTrackReference *tref)
398{
399 //
400 // Add track reference
401 //
402
eb05d549 403 if(fMC->fNTrackRefs >= 2*AliTRDgeometry::kNlayer){
1ee39b3a 404 SetCurved();
405 return;
406 }
407 // Make a copy for the object in order to avoid ownership problems
408 fMC->fTrackRefs[fMC->fNTrackRefs++] = new AliTrackReference(*tref);
409}
410
411//___________________________________________________
412AliTrackReference* AliTRDtrackInfo::GetTrackRef(Int_t idx) const
413{
414//
415// Returns a track reference
416//
f232df0d 417 if(!fMC) return NULL;
418 return (idx>=0 && idx < 12) ? fMC->fTrackRefs[idx] : NULL;
1ee39b3a 419}
420
421//___________________________________________________
61f6b45e 422AliTrackReference* AliTRDtrackInfo::GetTrackRef(const AliTRDseedV1* const tracklet) const
1ee39b3a 423{
424//
425// Returns a track reference
426//
f232df0d 427 if(!fMC) return NULL;
1ee39b3a 428 Double_t cw = AliTRDgeometry::CamHght() + AliTRDgeometry::CdrHght();
429 AliTrackReference * const* jtr = &(fMC->fTrackRefs[0]);
430 for(Int_t itr = 0; itr < fMC->fNTrackRefs; itr++, ++jtr){
431 if(!(*jtr)) break;
432 if(TMath::Abs(tracklet->GetX0() - (*jtr)->LocalX()) < cw) return (*jtr);
433 }
f232df0d 434 return NULL;
1ee39b3a 435}
436
437//___________________________________________________
438Int_t AliTRDtrackInfo::GetNumberOfClusters() const
439{
440 //
441 // Returns the number of clusters
442 //
443
444 Int_t n = 0;
445 if(!fTRDtrack) return 0;
446 if(fTRDtrack->GetNumberOfTracklets() == 0) return n;
f232df0d 447 AliTRDseedV1 *tracklet = NULL;
eb05d549 448 for(Int_t ip=0; ip<AliTRDgeometry::kNlayer; ip++){
1ee39b3a 449 if(!(tracklet = const_cast<AliTRDseedV1 *>(fTRDtrack->GetTracklet(ip)))) continue;
450 n+=tracklet->GetN();
451 }
452 return n;
453}
454
455
456//___________________________________________________
eb05d549 457void AliTRDtrackInfo::AliESDinfo::SetOuterParam(const AliExternalTrackParam *op)
1ee39b3a 458{
459 //
460 // Set outer track parameters
461 //
462
eb05d549 463 if(op){
464 if(fOP){
465 fOP->~AliExternalTrackParam();
466 // RS: Constructor from VTrack was used instead of Constructor from AliExternalTrackParam
467 new(fOP) AliExternalTrackParam(*op);
468 } else fOP = new AliExternalTrackParam(*op);
469 } else {
470 if(fOP) delete fOP; fOP = NULL;
471 }
1ee39b3a 472}
473
474//___________________________________________________
eb05d549 475void AliTRDtrackInfo::AliESDinfo::SetITSoutParam(const AliExternalTrackParam *op)
b9058d72 476{
477 //
478 // Set TPCout track parameters
479 //
480
eb05d549 481 if(op){
482 if(fITSout){
483 fITSout->~AliExternalTrackParam();
484 // RS: Constructor from VTrack was used instead of Constructor from AliExternalTrackParam
485 new(fITSout) AliExternalTrackParam(*op);
486 } else fITSout = new AliExternalTrackParam(*op);
487 } else {
488 if(fITSout) delete fITSout; fITSout = NULL;
489 }
b9058d72 490}
491
492//___________________________________________________
eb05d549 493void AliTRDtrackInfo::AliESDinfo::SetTPCoutParam(const AliExternalTrackParam *op)
35983729 494{
495 //
496 // Set TPCout track parameters
497 //
498
eb05d549 499 if(op){
500 if(fTPCout){
501 fTPCout->~AliExternalTrackParam();
502 // RS: Constructor from VTrack was used instead of Constructor from AliExternalTrackParam
503 new(fTPCout) AliExternalTrackParam(*op);
504 } else fTPCout = new AliExternalTrackParam(*op);
505 } else {
506 if(fTPCout) delete fTPCout; fTPCout = NULL;
507 }
35983729 508}
509
510//___________________________________________________
1ee39b3a 511Int_t AliTRDtrackInfo::GetNTracklets() const
512{
513 //
514 // Return the number of tracklets
515 //
516
e2927481 517 if(!fTRDtrack) return 0;
1ee39b3a 518 return fTRDtrack->GetNumberOfTracklets();
519}
520
521//___________________________________________________
522void AliTRDtrackInfo::SetSlices(Int_t n, Double32_t *s)
523{
524 //
525 // Set the slices
526 //
1ee39b3a 527 if(fESD.fTRDnSlices != n){
528 fESD.fTRDnSlices = 0;
529 delete [] fESD.fTRDslices;
f232df0d 530 fESD.fTRDslices = NULL;
1ee39b3a 531 }
532
533 if(!fESD.fTRDnSlices){
534 fESD.fTRDnSlices = n;
f232df0d 535 fESD.fTRDslices = new Double32_t[fESD.fTRDnSlices];
1ee39b3a 536 }
537
538 memcpy(fESD.fTRDslices, s, n*sizeof(Double32_t));
539}
540
541//___________________________________________________
0f2c4c4f 542Bool_t AliTRDtrackInfo::AliMCinfo::GetDirections(Float_t &x0, Float_t &y0, Float_t &z0, Float_t &dydx, Float_t &dzdx, Float_t &pt, Float_t &p, Float_t &eta, Float_t &phi, UChar_t &status) const
1ee39b3a 543{
544// Check for 2 track ref for the tracklet defined bythe radial position x0
545// The "status" is a bit map and gives a more informative output in case of failure:
546// - 0 : everything is OK
547// - BIT(0) : 0 track refs found
548// - BIT(1) : 1 track reference found
549// - BIT(2) : dx <= 0 between track references
550// - BIT(3) : dx > 0 && dx < 3.7 - tangent tracks
551
552 status = 0;
553 Double_t cw = AliTRDgeometry::CamHght() + AliTRDgeometry::CdrHght();
554 Int_t nFound = 0;
f232df0d 555 AliTrackReference *tr[2] = {NULL, NULL};
1ee39b3a 556 AliTrackReference * const* jtr = &fTrackRefs[0];
557 for(Int_t itr = 0; itr < fNTrackRefs; itr++, ++jtr){
558 if(!(*jtr)) break;
559/*
560 if(fDebugLevel>=5) printf("\t\tref[%2d] x[%6.3f]\n", itr, (*jtr)->LocalX());*/
561 if(TMath::Abs(x0 - (*jtr)->LocalX()) > cw) continue;
562 tr[nFound++] = (*jtr);
563 if(nFound == 2) break;
564 }
565 if(nFound < 2){
e02d9eb3 566 //AliDebug(1, Form("Missing track ref x0[%6.3f] nref[%d]", x0, nFound));
1ee39b3a 567 if(!nFound) SETBIT(status, 0);
568 else SETBIT(status, 1);
569 return kFALSE;
570 }
61f6b45e 571 pt=tr[1]->Pt();
0f2c4c4f 572 p =tr[1]->P();
61f6b45e 573 if(pt < 1.e-3) return kFALSE;
1ee39b3a 574
575 Double_t dx = tr[1]->LocalX() - tr[0]->LocalX();
576 if(dx <= 0.){
577 AliWarningGeneral("AliTRDtrackInfo::AliMCinfo::GetDirections()", Form("Track ref with wrong radial distances refX0[%6.3f] refX1[%6.3f]", tr[0]->LocalX(), tr[1]->LocalX()));
578 SETBIT(status, 2);
579 return kFALSE;
580 }
581 if(TMath::Abs(dx-AliTRDgeometry::CamHght()-AliTRDgeometry::CdrHght())>1.E-3) SETBIT(status, 3);
582
583 dydx = (tr[1]->LocalY() - tr[0]->LocalY()) / dx;
584 dzdx = (tr[1]->Z() - tr[0]->Z()) / dx;
585 //Float_t dx0 = tr[1]->LocalX() - x0;
586 y0 = tr[1]->LocalY()/* - dydx*dx0*/;
587 z0 = tr[1]->Z()/* - dzdx*dx0*/;
588 x0 = tr[1]->LocalX();
3ceb45ae 589 eta = -TMath::Log(TMath::Tan(0.5 * tr[1]->Theta()));
01ccc21a 590 phi = TMath::ATan2(tr[1]->Y(), tr[1]->X());
1ee39b3a 591 return kTRUE;
592}
593
594//___________________________________________________
3ceb45ae 595Bool_t AliTRDtrackInfo::AliMCinfo::PropagateKalman(
596 TVectorD *x, TVectorD *y, TVectorD *z,
597 TVectorD *dx, TVectorD *dy, TVectorD *dz,
598 TVectorD *pt, TVectorD *dpt, TVectorD *budget, TVectorD *c, Double_t mass) const
1ee39b3a 599{
600// Propagate Kalman from the first TRD track reference to
601// last one and save residuals in the y, z and pt.
602//
603// This is to calibrate the dEdx and MS corrections
604
3ceb45ae 605 if(!fNTrackRefs) return kFALSE;
1ee39b3a 606 for(Int_t itr=kNTrackRefs; itr--;){
3ceb45ae 607 (*x)[itr] = 0.;(*y)[itr] = 0.;(*z)[itr] = 0.;
4226db3e 608 (*dx)[itr] = -1.; (*dy)[itr] = 100.; (*dz)[itr] = 100.; (*dpt)[itr] = 100.;
1ee39b3a 609 }
1ee39b3a 610
611 // Initialize TRD track to the first track reference
4226db3e 612 AliTrackReference *tr(NULL);
613 Int_t itr(0); while(!(tr = fTrackRefs[itr])) itr++;
3ceb45ae 614 if(tr->Pt()<1.e-3) return kFALSE;
1ee39b3a 615
616 AliTRDtrackV1 tt;
617 Double_t xyz[3]={tr->X(),tr->Y(),tr->Z()};
618 Double_t pxyz[3]={tr->Px(),tr->Py(),tr->Pz()};
619 Double_t var[6] = {1.e-4, 1.e-4, 1.e-4, 1.e-4, 1.e-4, 1.e-4};
620 Double_t cov[21]={
621 var[0], 0., 0., 0., 0., 0.,
622 var[1], 0., 0., 0., 0.,
623 var[2], 0., 0., 0.,
624 var[3], 0., 0.,
625 var[4], 0.,
626 var[5]
627 };
628 TDatabasePDG db;
629 const TParticlePDG *pdg=db.GetParticle(fPDG);
630 if(!pdg){
631 AliWarningGeneral("AliTRDtrackInfo::AliMCinfo::PropagateKalman()", Form("PDG entry missing for code %d. References for track %d", fPDG, fNTrackRefs));
3ceb45ae 632 return kFALSE;
1ee39b3a 633 }
634 tt.Set(xyz, pxyz, cov, Short_t(pdg->Charge()));
3ceb45ae 635 if(mass<0){ // mass 0 use PDG
636 tt.SetMass(pdg->Mass());
637 } else { // use rec value
638 tt.SetMass(mass);
639 }
640
641// Double_t bg(tr->P()/pdg->Mass());
642// printf("\n\nNEW track PDG[%d] bg[%f] x[%f]\n", fPDG, bg, TMath::Log(bg));
4226db3e 643 Double_t x0(tr->LocalX());
644 const Double_t *cc(NULL);
645 for(Int_t ip=0; itr<fNTrackRefs; itr++){
646 if(!(tr = fTrackRefs[itr])) continue;
3ceb45ae 647// printf("ip[%d] det[%d]\n", ip, tr->DetectorId());
5066aa9a 648 if(!AliTRDtrackerV1::PropagateToX(tt, tr->LocalX(), fgKalmanStep)) continue;
1ee39b3a 649
650 //if(update) ...
3ceb45ae 651 tt.GetXYZ(xyz);
652 (*x)[ip] = xyz[0];
653 (*y)[ip] = xyz[1];
654 (*z)[ip] = xyz[2];
4226db3e 655 (*dx)[ip] = tt.GetX() - x0;
656 (*dy)[ip] = tt.GetY() - tr->LocalY();
657 (*dz)[ip] = tt.GetZ() - tr->Z();
5066aa9a 658 (*pt)[ip] = tr->Pt();
4226db3e 659 (*dpt)[ip] = tt.Pt()- tr->Pt();
3ceb45ae 660// printf("pt : kalman[%e] MC[%e]\n", tt.Pt(), tr->Pt());
661 (*budget)[ip] = tt.GetBudget(0);
1ee39b3a 662 cc = tt.GetCovariance();
3ceb45ae 663// printf("dx[%5.2f] sy[%f]\n", (*dx)[ip], TMath::Sqrt(cc[0]));
664 for(Int_t ic(0), jp(ip*15); ic<15; ic++, jp++) (*c)[jp]=cc[ic];
1ee39b3a 665 ip++;
666 }
3ceb45ae 667 return kTRUE;
1ee39b3a 668}