]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG1/TRD/info/AliTRDtrackInfo.cxx
major TRD reconstruction update
[u/mrichter/AliRoot.git] / PWG1 / 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"
33#include "AliExternalTrackParam.h"
34
35#include "AliTRDseedV1.h"
36#include "AliTRDtrackV1.h"
37#include "AliTRDgeometry.h"
38#include "AliTRDtrackerV1.h"
39
40#include "AliTRDtrackInfo.h"
41
42ClassImp(AliTRDtrackInfo)
43ClassImp(AliTRDtrackInfo::AliMCinfo)
44ClassImp(AliTRDtrackInfo::AliESDinfo)
5066aa9a 45Double_t AliTRDtrackInfo::AliMCinfo::fgKalmanStep = 2.;
1ee39b3a 46
47//___________________________________________________
48AliTRDtrackInfo::AliTRDtrackInfo():
49 TObject()
50 ,fNClusters(0)
f232df0d 51 ,fTRDtrack(NULL)
52 ,fMC(NULL)
1ee39b3a 53 ,fESD()
54{
55 //
56 // Default constructor
57 //
58}
59
60
61//___________________________________________________
62AliTRDtrackInfo::AliTRDtrackInfo(const AliTRDtrackInfo &trdInfo):
63 TObject((const TObject&)trdInfo)
64 ,fNClusters(trdInfo.fNClusters)
f232df0d 65 ,fTRDtrack(NULL)
a340fe39 66 ,fMC(NULL)
1ee39b3a 67 ,fESD(trdInfo.fESD)
68{
69 //
70 // copy Entries
71 //
72
73 if(trdInfo.fMC) fMC = new AliMCinfo(*trdInfo.fMC);
74
75 if(trdInfo.fTRDtrack){
76 fTRDtrack = new AliTRDtrackV1(*trdInfo.fTRDtrack);
77 if(trdInfo.fTRDtrack->IsOwner()) fTRDtrack->SetOwner();
78 }
79}
80
81//___________________________________________________
82AliTRDtrackInfo::AliMCinfo::AliMCinfo()
83 :fLabel(0)
84 ,fPDG(0)
85 ,fNTrackRefs(0)
86{
87 // Set 0-Pointers
88 memset(fTrackRefs, 0, sizeof(AliTrackReference *) * 12);
89}
90
91//___________________________________________________
92AliTRDtrackInfo::AliMCinfo::AliMCinfo(const AliMCinfo &mc)
93 :fLabel(mc.fLabel)
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)
135 ,fTRDpidQuality(0)
136 ,fTRDnSlices(0)
f232df0d 137 ,fTRDslices(NULL)
138 ,fOP(NULL)
1ee39b3a 139{
140 //
141 // Constructor
142 //
143
144 memset(fTRDr, 0, AliPID::kSPECIES*sizeof(Double32_t));
d80a6a00 145 memset(fTRDv0pid, 0, AliPID::kSPECIES*sizeof(Int_t));
1ee39b3a 146}
147
148//___________________________________________________
149AliTRDtrackInfo::AliESDinfo::AliESDinfo(const AliESDinfo &esd)
b9ddd472 150 :fHasV0(esd.fHasV0)
151 ,fId(esd.fId)
1ee39b3a 152 ,fStatus(esd.fStatus)
153 ,fKinkIndex(esd.fKinkIndex)
154 ,fTPCncls(esd.fTPCncls)
155 ,fTRDpidQuality(esd.fTRDpidQuality)
156 ,fTRDnSlices(esd.fTRDnSlices)
f232df0d 157 ,fTRDslices(NULL)
158 ,fOP(NULL)
1ee39b3a 159{
160 //
161 // Constructor
162 //
163
164 memcpy(fTRDr, esd.fTRDr, AliPID::kSPECIES*sizeof(Double32_t));
d80a6a00 165 memcpy(fTRDv0pid, esd.fTRDv0pid, AliPID::kSPECIES*sizeof(Int_t));
1ee39b3a 166
167 if(fTRDnSlices){
168 fTRDslices = new Double32_t[fTRDnSlices];
169 memcpy(fTRDslices, esd.fTRDslices, fTRDnSlices*sizeof(Double32_t));
170 }
171 if(esd.fOP) fOP = new AliExternalTrackParam(*esd.fOP);
172}
173
174
175//___________________________________________________
176AliTRDtrackInfo::~AliTRDtrackInfo()
177{
178 //
179 // Destructor
180 //
181
182 if(fMC) delete fMC;
183 if(fTRDtrack) delete fTRDtrack;
184}
185
186//___________________________________________________
187AliTRDtrackInfo::AliMCinfo::~AliMCinfo()
188{
189 //
190 // Destructor
191 //
192
193 fNTrackRefs = 0;
194 for(Int_t ien = 0; ien < 12; ien++){
195 if(fTrackRefs[ien]) delete fTrackRefs[ien];
f232df0d 196 fTrackRefs[ien] = NULL;
1ee39b3a 197 }
198}
199
200//___________________________________________________
201AliTRDtrackInfo::AliESDinfo::~AliESDinfo()
202{
203 //
204 // Destructor
205 //
206
207 if(fTRDnSlices){
208 delete [] fTRDslices;
f232df0d 209 fTRDslices = NULL;
1ee39b3a 210 fTRDnSlices = 0;
211 }
f232df0d 212 if(fOP) delete fOP; fOP = NULL;
1ee39b3a 213}
214
8778aca3 215//___________________________________________________
216void AliTRDtrackInfo::AliESDinfo::Delete(const Option_t *){
217 //
218 // Delete Pointer members
219 //
220 if(fTRDnSlices){
221 delete [] fTRDslices;
222 fTRDslices = NULL;
223 fTRDnSlices = 0;
224 }
225 if(fOP) delete fOP; fOP = NULL;
226}
227
1ee39b3a 228
229//___________________________________________________
230AliTRDtrackInfo& AliTRDtrackInfo::operator=(const AliTRDtrackInfo &trdInfo)
231{
232 //
233 // = Operator
234 //
235
236 fNClusters = trdInfo.fNClusters;
237 fESD = trdInfo.fESD;
238
239 if(trdInfo.fMC){
2c33fb46 240 if(!fMC) fMC = new AliMCinfo(*trdInfo.fMC);
241 else{
242 fMC->~AliMCinfo();
1ee39b3a 243 new(fMC) AliMCinfo(*trdInfo.fMC);
2c33fb46 244 }
1ee39b3a 245 }
246
247 if(trdInfo.fTRDtrack){
2c33fb46 248 if(!fTRDtrack) fTRDtrack = new AliTRDtrackV1(*trdInfo.fTRDtrack);
249 else{
250 fTRDtrack->~AliTRDtrackV1();
1ee39b3a 251 new(fTRDtrack) AliTRDtrackV1(*trdInfo.fTRDtrack);
2c33fb46 252 }
1ee39b3a 253 if(trdInfo.fTRDtrack->IsOwner()) fTRDtrack->SetOwner();
254 }
255
256 return *this;
257}
258
259//___________________________________________________
260AliTRDtrackInfo::AliMCinfo& AliTRDtrackInfo::AliMCinfo::operator=(const AliMCinfo &mc)
261{
262 //
263 // Assignment operator
264 //
265
266 fLabel = mc.fLabel;
267 fPDG = mc.fPDG;
268 fNTrackRefs = mc.fNTrackRefs;
269
270 AliTrackReference **itr = &fTrackRefs[0];
271 AliTrackReference* const *jtr = &mc.fTrackRefs[0];
272 for(Int_t ien = 0; ien < 12; ien++, itr++, jtr++){
273 if((*jtr)){
274 if(!(*itr)) (*itr) = new AliTrackReference(*(*jtr));
2c33fb46 275 else{
276 (*itr)->~AliTrackReference();
277 new(&(*itr)) AliTrackReference(*(*jtr));
278 }
f232df0d 279 } else (*itr) = NULL;
1ee39b3a 280 }
281 return *this;
282}
283
284//___________________________________________________
285AliTRDtrackInfo::AliESDinfo& AliTRDtrackInfo::AliESDinfo::operator=(const AliESDinfo &esd)
286{
287 //
288 // Assignment operator
289 //
290
291 fId = esd.fId;
292 fStatus = esd.fStatus;
293 fTRDpidQuality = esd.fTRDpidQuality;
294 fTRDnSlices = esd.fTRDnSlices;
f232df0d 295 fTRDslices = NULL;
1ee39b3a 296
297 memcpy(fTRDr, esd.fTRDr, AliPID::kSPECIES*sizeof(Double32_t));
298
299 if(fTRDnSlices){
300 fTRDslices = new Double32_t[fTRDnSlices];
301 memcpy(fTRDslices, esd.fTRDslices, fTRDnSlices*sizeof(Double32_t));
302 }
303 if(esd.fOP){
2c33fb46 304 if(fOP){
305 fOP->~AliExternalTrackParam();
306 new(fOP) AliExternalTrackParam(esd.fOP);
307 } else fOP = new AliExternalTrackParam(esd.fOP);
f232df0d 308 } else fOP = NULL;
1ee39b3a 309
310 return *this;
311}
312
313//___________________________________________________
314void AliTRDtrackInfo::Delete(const Option_t *)
315{
316 //
317 // Delete
318 //
319
2c33fb46 320 AliDebug(2, Form("track[%p] mc[%p]", (void*)fTRDtrack, (void*)fMC));
1ee39b3a 321 fNClusters = 0;
f232df0d 322 if(fMC) delete fMC; fMC = NULL;
8778aca3 323 fESD.Delete(NULL);
f232df0d 324 if(fTRDtrack) delete fTRDtrack; fTRDtrack = NULL;
1ee39b3a 325}
326
327//___________________________________________________
328void AliTRDtrackInfo::SetTrack(const AliTRDtrackV1 *track)
329{
330 //
331 // Set the TRD track
332 //
333
2c33fb46 334 if(!fTRDtrack) fTRDtrack = new AliTRDtrackV1(*track);
335 else{
336 fTRDtrack->~AliTRDtrackV1();
337 new(fTRDtrack) AliTRDtrackV1(*track);
338 }
1ee39b3a 339 fTRDtrack->SetOwner();
340 // Make a copy for the object in order to avoid ownership problems
341}
342
343//___________________________________________________
344void AliTRDtrackInfo::AddTrackRef(const AliTrackReference *tref)
345{
346 //
347 // Add track reference
348 //
349
1ee39b3a 350 if(fMC->fNTrackRefs >= 12){
351 SetCurved();
352 return;
353 }
354 // Make a copy for the object in order to avoid ownership problems
355 fMC->fTrackRefs[fMC->fNTrackRefs++] = new AliTrackReference(*tref);
356}
357
358//___________________________________________________
359AliTrackReference* AliTRDtrackInfo::GetTrackRef(Int_t idx) const
360{
361//
362// Returns a track reference
363//
f232df0d 364 if(!fMC) return NULL;
365 return (idx>=0 && idx < 12) ? fMC->fTrackRefs[idx] : NULL;
1ee39b3a 366}
367
368//___________________________________________________
61f6b45e 369AliTrackReference* AliTRDtrackInfo::GetTrackRef(const AliTRDseedV1* const tracklet) const
1ee39b3a 370{
371//
372// Returns a track reference
373//
f232df0d 374 if(!fMC) return NULL;
1ee39b3a 375 Double_t cw = AliTRDgeometry::CamHght() + AliTRDgeometry::CdrHght();
376 AliTrackReference * const* jtr = &(fMC->fTrackRefs[0]);
377 for(Int_t itr = 0; itr < fMC->fNTrackRefs; itr++, ++jtr){
378 if(!(*jtr)) break;
379 if(TMath::Abs(tracklet->GetX0() - (*jtr)->LocalX()) < cw) return (*jtr);
380 }
f232df0d 381 return NULL;
1ee39b3a 382}
383
384//___________________________________________________
385Int_t AliTRDtrackInfo::GetNumberOfClusters() const
386{
387 //
388 // Returns the number of clusters
389 //
390
391 Int_t n = 0;
392 if(!fTRDtrack) return 0;
393 if(fTRDtrack->GetNumberOfTracklets() == 0) return n;
f232df0d 394 AliTRDseedV1 *tracklet = NULL;
1ee39b3a 395 for(Int_t ip=0; ip<6; ip++){
396 if(!(tracklet = const_cast<AliTRDseedV1 *>(fTRDtrack->GetTracklet(ip)))) continue;
397 n+=tracklet->GetN();
398 }
399 return n;
400}
401
402
403//___________________________________________________
404void AliTRDtrackInfo::SetOuterParam(const AliExternalTrackParam *op)
405{
406 //
407 // Set outer track parameters
408 //
409
410 if(!op) return;
2c33fb46 411 if(fESD.fOP){
412 fESD.fOP->~AliExternalTrackParam();
413 new(fESD.fOP) AliExternalTrackParam(*op);
414 } else fESD.fOP = new AliExternalTrackParam(*op);
1ee39b3a 415}
416
417//___________________________________________________
418Int_t AliTRDtrackInfo::GetNTracklets() const
419{
420 //
421 // Return the number of tracklets
422 //
423
e2927481 424 if(!fTRDtrack) return 0;
1ee39b3a 425 return fTRDtrack->GetNumberOfTracklets();
426}
427
428//___________________________________________________
429void AliTRDtrackInfo::SetSlices(Int_t n, Double32_t *s)
430{
431 //
432 // Set the slices
433 //
1ee39b3a 434 if(fESD.fTRDnSlices != n){
435 fESD.fTRDnSlices = 0;
436 delete [] fESD.fTRDslices;
f232df0d 437 fESD.fTRDslices = NULL;
1ee39b3a 438 }
439
440 if(!fESD.fTRDnSlices){
441 fESD.fTRDnSlices = n;
f232df0d 442 fESD.fTRDslices = new Double32_t[fESD.fTRDnSlices];
1ee39b3a 443 }
444
445 memcpy(fESD.fTRDslices, s, n*sizeof(Double32_t));
446}
447
448//___________________________________________________
449Bool_t AliTRDtrackInfo::AliMCinfo::GetDirections(Float_t &x0, Float_t &y0, Float_t &z0, Float_t &dydx, Float_t &dzdx, Float_t &pt, UChar_t &status) const
450{
451// Check for 2 track ref for the tracklet defined bythe radial position x0
452// The "status" is a bit map and gives a more informative output in case of failure:
453// - 0 : everything is OK
454// - BIT(0) : 0 track refs found
455// - BIT(1) : 1 track reference found
456// - BIT(2) : dx <= 0 between track references
457// - BIT(3) : dx > 0 && dx < 3.7 - tangent tracks
458
459 status = 0;
460 Double_t cw = AliTRDgeometry::CamHght() + AliTRDgeometry::CdrHght();
461 Int_t nFound = 0;
f232df0d 462 AliTrackReference *tr[2] = {NULL, NULL};
1ee39b3a 463 AliTrackReference * const* jtr = &fTrackRefs[0];
464 for(Int_t itr = 0; itr < fNTrackRefs; itr++, ++jtr){
465 if(!(*jtr)) break;
466/*
467 if(fDebugLevel>=5) printf("\t\tref[%2d] x[%6.3f]\n", itr, (*jtr)->LocalX());*/
468 if(TMath::Abs(x0 - (*jtr)->LocalX()) > cw) continue;
469 tr[nFound++] = (*jtr);
470 if(nFound == 2) break;
471 }
472 if(nFound < 2){
f232df0d 473 AliWarningGeneral("AliTRDtrackInfo::AliMCinfo::GetDirections()", Form("Missing track ref x0[%6.3f] nref[%d]", x0, nFound));
1ee39b3a 474 if(!nFound) SETBIT(status, 0);
475 else SETBIT(status, 1);
476 return kFALSE;
477 }
61f6b45e 478 pt=tr[1]->Pt();
479 if(pt < 1.e-3) return kFALSE;
1ee39b3a 480
481 Double_t dx = tr[1]->LocalX() - tr[0]->LocalX();
482 if(dx <= 0.){
483 AliWarningGeneral("AliTRDtrackInfo::AliMCinfo::GetDirections()", Form("Track ref with wrong radial distances refX0[%6.3f] refX1[%6.3f]", tr[0]->LocalX(), tr[1]->LocalX()));
484 SETBIT(status, 2);
485 return kFALSE;
486 }
487 if(TMath::Abs(dx-AliTRDgeometry::CamHght()-AliTRDgeometry::CdrHght())>1.E-3) SETBIT(status, 3);
488
489 dydx = (tr[1]->LocalY() - tr[0]->LocalY()) / dx;
490 dzdx = (tr[1]->Z() - tr[0]->Z()) / dx;
491 //Float_t dx0 = tr[1]->LocalX() - x0;
492 y0 = tr[1]->LocalY()/* - dydx*dx0*/;
493 z0 = tr[1]->Z()/* - dzdx*dx0*/;
494 x0 = tr[1]->LocalX();
495 return kTRUE;
496}
497
498//___________________________________________________
5066aa9a 499void AliTRDtrackInfo::AliMCinfo::PropagateKalman(TVectorD *dx, TVectorD *dy, TVectorD *dz, TVectorD *pt, TVectorD *dpt, TVectorD *c) const
1ee39b3a 500{
501// Propagate Kalman from the first TRD track reference to
502// last one and save residuals in the y, z and pt.
503//
504// This is to calibrate the dEdx and MS corrections
505
4226db3e 506 if(!fNTrackRefs) return;
1ee39b3a 507 for(Int_t itr=kNTrackRefs; itr--;){
4226db3e 508 (*dx)[itr] = -1.; (*dy)[itr] = 100.; (*dz)[itr] = 100.; (*dpt)[itr] = 100.;
1ee39b3a 509 }
1ee39b3a 510
511 // Initialize TRD track to the first track reference
4226db3e 512 AliTrackReference *tr(NULL);
513 Int_t itr(0); while(!(tr = fTrackRefs[itr])) itr++;
1ee39b3a 514 if(tr->Pt()<1.e-3) return;
515
516 AliTRDtrackV1 tt;
517 Double_t xyz[3]={tr->X(),tr->Y(),tr->Z()};
518 Double_t pxyz[3]={tr->Px(),tr->Py(),tr->Pz()};
519 Double_t var[6] = {1.e-4, 1.e-4, 1.e-4, 1.e-4, 1.e-4, 1.e-4};
520 Double_t cov[21]={
521 var[0], 0., 0., 0., 0., 0.,
522 var[1], 0., 0., 0., 0.,
523 var[2], 0., 0., 0.,
524 var[3], 0., 0.,
525 var[4], 0.,
526 var[5]
527 };
528 TDatabasePDG db;
529 const TParticlePDG *pdg=db.GetParticle(fPDG);
530 if(!pdg){
531 AliWarningGeneral("AliTRDtrackInfo::AliMCinfo::PropagateKalman()", Form("PDG entry missing for code %d. References for track %d", fPDG, fNTrackRefs));
532 return;
533 }
534 tt.Set(xyz, pxyz, cov, Short_t(pdg->Charge()));
535 tt.SetMass(pdg->Mass());
536
4226db3e 537 Double_t x0(tr->LocalX());
538 const Double_t *cc(NULL);
539 for(Int_t ip=0; itr<fNTrackRefs; itr++){
540 if(!(tr = fTrackRefs[itr])) continue;
5066aa9a 541 if(!AliTRDtrackerV1::PropagateToX(tt, tr->LocalX(), fgKalmanStep)) continue;
1ee39b3a 542
543 //if(update) ...
4226db3e 544 (*dx)[ip] = tt.GetX() - x0;
545 (*dy)[ip] = tt.GetY() - tr->LocalY();
546 (*dz)[ip] = tt.GetZ() - tr->Z();
5066aa9a 547 (*pt)[ip] = tr->Pt();
4226db3e 548 (*dpt)[ip] = tt.Pt()- tr->Pt();
1ee39b3a 549 cc = tt.GetCovariance();
4226db3e 550 c->Use(ip*15, (ip+1)*15, cc);
1ee39b3a 551 ip++;
552 }
553}