PWG1/TRD/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 }