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