PWG1/TPC/AliESDRecInfo.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
  17 ///////////////////////////////////////////////////////////////////////////////
  18 //                                                                           //
  19 //  Time Projection Chamber                                                  //
  20 //  Comparison macro for reconstructed tracks                                  //
  21 //  responsible:
  22 //  marian.ivanov@cern.ch                                                    //
  23 //
  24 //
  25
  26
  27
  28
  29
  30 //ROOT includes
  31 #include "Rtypes.h"
  32 //
  33 //ALIROOT includes
  34 //
  35 #include "AliESDtrack.h"
  36 #include "AliTracker.h"
  37 #include "AliTPCParam.h"
  38 #include "AliTrackReference.h"
  39 #include "AliTPCParamSR.h"
  40 #include "AliESDfriend.h"
  41 #include "AliESDtrack.h"
  42 #include "AliTPCseed.h"
  43 #include "AliITStrackMI.h"
  44 #include "AliTRDtrackV1.h"
  45 #include "AliMCInfo.h"
  46 #include "AliESDRecInfo.h"
  47
  48
  49
  50 ClassImp(AliESDRecInfo)
  51
  52
  53
  54
  55 AliTPCParam * GetTPCParam(){
  56   AliTPCParamSR * par = new AliTPCParamSR;
  57   par->Update();
  58   return par;
  59 }
  60
  61
  62
  63
  64 AliESDRecInfo::AliESDRecInfo():
  65   fITSOn(0),           // ITS refitted inward
  66   fTRDOn(0),           // ITS refitted inward
  67   fDeltaP(0),          //delta of momenta
  68   fSign(0),           // sign
  69   fReconstructed(0),         //flag if track was reconstructed
  70   fFake(0),             // fake track
  71   fMultiple(0),         // number of reconstructions
  72   fTPCOn(0),           // TPC refitted inward
  73   fBestTOFmatch(0),        //best matching between times
  74   fESDtrack(0),        // esd track
  75   fTrackF(0),      // friend track
  76   fTPCtrack(0),        // tpc track
  77   fITStrack(0),        // its track
  78   fTRDtrack(0),        // trd track
  79   fTracks(0)           // array of tracks with the same label
  80 {
  81   //
  82   //  default constructor
  83   //
  84 }
  85
  86
  87 AliESDRecInfo::AliESDRecInfo(const AliESDRecInfo& recinfo):
  88   TObject(),
  89   fITSOn(0),           // ITS refitted inward
  90   fTRDOn(0),           // ITS refitted inward
  91   fDeltaP(0),          //delta of momenta
  92   fSign(0),           // sign
  93   fReconstructed(0),         //flag if track was reconstructed
  94   fFake(0),             // fake track
  95   fMultiple(0),         // number of reconstructions
  96   fTPCOn(0),           // TPC refitted inward
  97   fBestTOFmatch(0),        //best matching between times
  98   fESDtrack(0),        // esd track
  99   fTrackF(0),      // friend track
 100   fTPCtrack(0),        // tpc track
 101   fITStrack(0),        // its track
 102   fTRDtrack(0),        // trd track
 103   fTracks(0)           // array of tracks with the same label
 104 {
 105   //
 106   //
 107   //
 108   memcpy(this,&recinfo, sizeof(recinfo));
 109   fESDtrack=0; fTrackF=0; fTPCtrack=0;fITStrack=0;fTRDtrack=0;
 110   fTracks=0;
 111   SetESDtrack(recinfo.GetESDtrack());
 112 }
 113
 114
 115 AliESDRecInfo& AliESDRecInfo::operator=(const AliESDRecInfo& info) {
 116   //
 117   // Assignment operator
 118   //
 119   this->~AliESDRecInfo();
 120   new (this) AliESDRecInfo(info);
 121   return *this;
 122 }
 123
 124
 125
 126 AliESDRecInfo::~AliESDRecInfo()
 127
 128 {
 129   //
 130   //  destructor
 131   //
 132   if (fESDtrack) { delete fESDtrack; fESDtrack=0;}
 133   if (fTrackF)   {   fTrackF=0;}
 134   if (fTPCtrack) { delete fTPCtrack; fTPCtrack=0;}
 135   if (fITStrack) { delete fITStrack; fITStrack=0;}
 136   if (fTRDtrack) { delete fTRDtrack; fTRDtrack=0;}
 137   if (fTracks) {
 138     delete  fTracks;  fTracks=0;
 139   }
 140 }
 141
 142
 143
 144 void AliESDRecInfo::Reset()
 145 {
 146   //
 147   // reset info
 148   //
 149   fMultiple =0;
 150   fFake     =0;
 151   fReconstructed=0;
 152   if (fESDtrack) { delete fESDtrack; fESDtrack=0;}
 153   if (fTrackF)   { fTrackF=0;}
 154   if (fTPCtrack) { delete fTPCtrack; fTPCtrack=0;}
 155   if (fITStrack) { delete fITStrack; fITStrack=0;}
 156   if (fTRDtrack) { delete fTRDtrack; fTRDtrack=0;}
 157   if (fTracks) { delete  fTracks;  fTracks=0;}
 158 }
 159
 160 void AliESDRecInfo::SetESDtrack(const AliESDtrack *track){
 161   //
 162   //
 163   //
 164
 165   if (fESDtrack) delete fESDtrack;
 166   fESDtrack = (AliESDtrack*)track->Clone();
 167   //AliESDfriendTrack *friendTrack=fESDfriend->GetTrack(fESDtrack->GetID());
 168
 169
 170   if (track->GetFriendTrack()){
 171     fTrackF = (AliESDfriendTrack*)track->GetFriendTrack();
 172     Int_t icalib=0;
 173     TObject *cobject=0;
 174     //
 175     while (fTrackF->GetCalibObject(icalib)){
 176       cobject=fTrackF->GetCalibObject(icalib);
 177       if (dynamic_cast<AliTPCseed*>(cobject)){
 178         if (fTPCtrack) delete fTPCtrack;
 179         fTPCtrack = (AliTPCseed*)(dynamic_cast<AliTPCseed*>(cobject))->Clone();
 180       } else if (dynamic_cast<AliTRDtrackV1*>(cobject)){
 181         if (fTRDtrack) delete fTRDtrack;
 182         fTRDtrack = (AliTRDtrackV1*)(dynamic_cast<AliTRDtrackV1*>(cobject))->Clone();
 183       }
 184       icalib++;
 185     }
 186   }
 187   if (!fTPCtrack) fTPCtrack = new AliTPCseed; // add dummy track
 188 }
 189
 190
 191
 192 void AliESDRecInfo::AddESDtrack(const AliESDtrack *track, AliMCInfo* info){
 193   //
 194   // Add ESDtrack
 195   //
 196   AliESDtrack *nctrack = (AliESDtrack*) track;
 197   fMultiple++;
 198   if (!fESDtrack) {
 199     SetESDtrack(track);
 200     Update(info,0,kTRUE);
 201     return;
 202   }
 203   if (!fTracks) fTracks = new TClonesArray("AliESDtrack",10);
 204   Int_t ntracks = fTracks->GetEntriesFast();
 205   new ((*fTracks)[ntracks]) AliESDtrack(*track);
 206   if  (nctrack->GetKinkIndex(0)>0) return;
 207   //
 208   //
 209   //
 210   if (!nctrack->IsOn(AliESDtrack::kTPCrefit) && fStatus[1]==2) return;
 211   if (!nctrack->IsOn(AliESDtrack::kITSin) && fStatus[0]>0) return;
 212   if ( nctrack->GetParameter()[4]*info->GetCharge()<0) return; //returning track
 213
 214   Float_t dtheta = TMath::ATan(nctrack->GetTgl())-info->GetParticle().Theta()-TMath::Pi()-2;
 215   if (TMath::Abs(dtheta)>0.1) return;
 216
 217   SetESDtrack(track);
 218   Update(info,0,kTRUE);
 219 }
 220
 221
 222
 223 void  AliESDRecInfo::UpdatePoints(AliESDtrack*track)
 224 {
 225   //
 226   //
 227   Int_t iclusters[200];
 228   Float_t density[160];
 229   for (Int_t i=0;i<160;i++) density[i]=-1.;
 230   fTPCPoints[0]= 160;
 231   fTPCPoints[1] = -1;
 232   //
 233   if (fTPCPoints[0]<fTPCPoints[1]) return;
 234     /*Int_t nclusters=*/track->GetTPCclusters(iclusters);
 235
 236   Int_t ngood=0;
 237   Int_t undeff=0;
 238   Int_t nall =0;
 239   Int_t range=20;
 240   for (Int_t i=0;i<160;i++){
 241     Int_t last = i-range;
 242     if (nall<range) nall++;
 243     if (last>=0){
 244       if (iclusters[last]>0&& (iclusters[last]&0x8000)==0) ngood--;
 245       if (iclusters[last]==-1) undeff--;
 246     }
 247     if (iclusters[i]>0&& (iclusters[i]&0x8000)==0)   ngood++;
 248     if (iclusters[i]==-1) undeff++;
 249     if (nall==range &&undeff<range/2) density[i-range/2] = Float_t(ngood)/Float_t(nall-undeff);
 250   }
 251   Float_t maxdens=0;
 252   Int_t indexmax =0;
 253   for (Int_t i=0;i<160;i++){
 254     if (density[i]<0) continue;
 255     if (density[i]>maxdens){
 256       maxdens=density[i];
 257       indexmax=i;
 258     }
 259   }
 260   //
 261   //max dens point
 262   fTPCPoints[3] = maxdens;
 263   fTPCPoints[1] = indexmax;
 264   //
 265   // last point
 266   for (Int_t i=indexmax;i<160;i++){
 267     if (density[i]<0) continue;
 268     if (density[i]<maxdens/2.) {
 269       break;
 270     }
 271     fTPCPoints[2]=i;
 272   }
 273   //
 274   // first point
 275   for (Int_t i=indexmax;i>0;i--){
 276     if (density[i]<0) continue;
 277     if (density[i]<maxdens/2.) {
 278       break;
 279     }
 280     fTPCPoints[0]=i;
 281   }
 282   //
 283   // Density at the last 30 padrows
 284   //
 285   //
 286   nall  = 0;
 287   ngood = 0;
 288   for (Int_t i=159;i>0;i--){
 289     if (iclusters[i]==-1) continue; //dead zone
 290     nall++;
 291     if (iclusters[i]>0)   ngood++;
 292     if (nall>20) break;
 293   }
 294   fTPCPoints[4] = Float_t(ngood)/Float_t(nall);
 295   //
 296   if ((track->GetStatus()&AliESDtrack::kITSrefit)>0) fTPCPoints[0]=-1;
 297
 298
 299 }
 300
 301 //
 302 //
 303 void AliESDRecInfo::Update(AliMCInfo* info,AliTPCParam * /*par*/, Bool_t reconstructed)
 304 {
 305   //
 306   //calculates derived variables
 307   //
 308   UpdatePoints(fESDtrack);
 309   UpdateStatus(info,reconstructed);
 310   UpdateITS(info);
 311   UpdateTPC(info);
 312   UpdateTOF(info);
 313 }
 314
 315
 316 void  AliESDRecInfo::UpdateStatus(AliMCInfo* info, Bool_t reconstructed){
 317   //
 318   // Interpret bit mask flags
 319   //
 320   for (Int_t i=0;i<4;i++) fStatus[i] =0;
 321   fReconstructed = kFALSE;
 322   fTPCOn = kFALSE;
 323   fITSOn = kFALSE;
 324   fTRDOn = kFALSE;
 325   if (reconstructed==kFALSE) return;
 326
 327   fLabels[0] = info->fLabel;
 328   fLabels[1] = info->fPrimPart;
 329   fReconstructed = kTRUE;
 330   fTPCOn = ((fESDtrack->GetStatus()&AliESDtrack::kTPCrefit)>0) ? kTRUE : kFALSE;
 331   fITSOn = ((fESDtrack->GetStatus()&AliESDtrack::kITSrefit)>0) ? kTRUE : kFALSE;
 332   fTRDOn = ((fESDtrack->GetStatus()&AliESDtrack::kTRDrefit)>0) ? kTRUE : kFALSE;
 333   //
 334   //
 335   if ((fESDtrack->GetStatus()&AliESDtrack::kTPCrefit)>0){
 336     fStatus[1] =3;
 337   }
 338   else{
 339     if ((fESDtrack->GetStatus()&AliESDtrack::kTPCout)>0){
 340       fStatus[1] =2;
 341     }
 342     else{
 343       if ((fESDtrack->GetStatus()&AliESDtrack::kTPCin)>0)
 344         fStatus[1]=1;
 345     }
 346   }
 347   //
 348   if ((fESDtrack->GetStatus()&AliESDtrack::kITSout)>0){
 349     fStatus[0] =2;
 350   }
 351   else{
 352     if ((fESDtrack->GetStatus()&AliESDtrack::kITSrefit)>0){
 353       fStatus[0] =1;
 354     }
 355     else{
 356       fStatus[0]=0;
 357     }
 358   }
 359   //
 360   //
 361   if ((fESDtrack->GetStatus()&AliESDtrack::kTRDrefit)>0){
 362     fStatus[2] =2;
 363   }
 364   else{
 365     if ((fESDtrack->GetStatus()&AliESDtrack::kTRDout)>0){
 366       fStatus[2] =1;
 367     }
 368   }
 369   if ((fESDtrack->GetStatus()&AliESDtrack::kTRDStop)>0){
 370     fStatus[2] =10;
 371   }
 372 }
 373
 374 void AliESDRecInfo::UpdateTOF(AliMCInfo* info){
 375   //
 376   // Update TOF related comparison information
 377   //
 378   fBestTOFmatch=1000;
 379   if (((fESDtrack->GetStatus()&AliESDtrack::kTOFout)>0)){
 380     //
 381     // best tof match
 382     Double_t times[5];
 383     fESDtrack->GetIntegratedTimes(times);
 384     for (Int_t i=0;i<5;i++){
 385       if ( TMath::Abs(fESDtrack->GetTOFsignal()-times[i]) <TMath::Abs(fBestTOFmatch) ){
 386         fBestTOFmatch = fESDtrack->GetTOFsignal()-times[i];
 387       }
 388     }
 389     Int_t toflabel[3];
 390     fESDtrack->GetTOFLabel(toflabel);
 391     Bool_t toffake=kTRUE;
 392     Bool_t tofdaughter=kFALSE;
 393     for (Int_t i=0;i<3;i++){
 394       if (toflabel[i]<0) continue;
 395       if (toflabel[i]== TMath::Abs(fESDtrack->GetLabel()))  toffake=kFALSE;
 396       if (toflabel[i]==info->fParticle.GetDaughter(0) || (toflabel[i]==info->fParticle.GetDaughter(1))) tofdaughter=kTRUE;  // decay product of original particle
 397       fStatus[3]=1;
 398     }
 399     if (toffake) fStatus[3] =3;       //total fake
 400     if (tofdaughter) fStatus[3]=2;    //fake because of decay
 401   }else{
 402     fStatus[3]=0;
 403   }
 404 }
 405
 406
 407
 408 void AliESDRecInfo::UpdateITS(AliMCInfo* info){
 409   //
 410   // Update ITS related comparison information
 411   //
 412   fITSinP0[0]=info->fParticle.Px();
 413   fITSinP0[1]=info->fParticle.Py();
 414   fITSinP0[2]=info->fParticle.Pz();
 415   fITSinP0[3]=info->fParticle.Pt();
 416   //
 417   fITSinR0[0]=info->fParticle.Vx();
 418   fITSinR0[1]=info->fParticle.Vy();
 419   fITSinR0[2]=info->fParticle.Vz();
 420   fITSinR0[3] = TMath::Sqrt(fITSinR0[0]*fITSinR0[0]+fITSinR0[1]*fITSinR0[1]);
 421   fITSinR0[4] = TMath::ATan2(fITSinR0[1],fITSinR0[0]);
 422   //
 423   //
 424   if (fITSinP0[3]>0.0000001){
 425     fITSAngle0[0] = TMath::ATan2(fITSinP0[1],fITSinP0[0]);
 426     fITSAngle0[1] = TMath::ATan(fITSinP0[2]/fITSinP0[3]);
 427   }
 428   if (fITSOn){
 429     // ITS
 430     Double_t param[5],x;
 431     fESDtrack->GetExternalParameters(x,param);
 432     //    fESDtrack->GetConstrainedExternalParameters(x,param);
 433     Double_t cov[15];
 434     fESDtrack->GetExternalCovariance(cov);
 435     //fESDtrack->GetConstrainedExternalCovariance(cov);
 436     if (TMath::Abs(param[4])<0.0000000001) return;
 437
 438     fESDtrack->GetXYZ(fITSinR1);
 439     fESDtrack->GetPxPyPz(fITSinP1);
 440     fITSinP1[3] = TMath::Sqrt(fITSinP1[0]*fITSinP1[0]+fITSinP1[1]*fITSinP1[1]);
 441     //
 442     fITSinR1[3] = TMath::Sqrt(fITSinR1[0]*fITSinR1[0]+fITSinR1[1]*fITSinR1[1]);
 443     fITSinR1[4] = TMath::ATan2(fITSinR1[1],fITSinR1[0]);
 444     //
 445     //
 446     if (fITSinP1[3]>0.0000001){
 447       fITSAngle1[0] = TMath::ATan2(fITSinP1[1],fITSinP1[0]);
 448       fITSAngle1[1] = TMath::ATan(fITSinP1[2]/fITSinP1[3]);
 449     }
 450     //
 451     //
 452     fITSDelta[0] = (fITSinR0[4]-fITSinR1[4])*fITSinR1[3];  //delta rfi
 453     fITSPools[0] = fITSDelta[0]/TMath::Sqrt(cov[0]);
 454     fITSDelta[1] = (fITSinR0[2]-fITSinR1[2]);              //delta z
 455     fITSPools[1] = fITSDelta[1]/TMath::Sqrt(cov[2]);
 456     fITSDelta[2] = (fITSAngle0[0]-fITSAngle1[0]);
 457     fITSPools[2] = fITSDelta[2]/TMath::Sqrt(cov[5]);
 458     fITSDelta[3] = (TMath::Tan(fITSAngle0[1])-TMath::Tan(fITSAngle1[1]));
 459     fITSPools[3] = fITSDelta[3]/TMath::Sqrt(cov[9]);
 460     fITSDelta[4] = (fITSinP0[3]-fITSinP1[3]);
 461     Double_t sign = (param[4]>0) ? 1:-1;
 462     fSign = sign;
 463     fITSPools[4] = sign*(1./fITSinP0[3]-1./fITSinP1[3])/TMath::Sqrt(cov[14]);
 464   }
 465 }
 466
 467 void AliESDRecInfo::UpdateTPC(AliMCInfo* info){
 468   //
 469   //
 470   // Update TPC related information
 471   //
 472   AliTrackReference * ref = &(info->fTrackRef);
 473   fTPCinR0[0] = info->fTrackRef.X();
 474   fTPCinR0[1] = info->fTrackRef.Y();
 475   fTPCinR0[2] = info->fTrackRef.Z();
 476   fTPCinR0[3] = TMath::Sqrt(fTPCinR0[0]*fTPCinR0[0]+fTPCinR0[1]*fTPCinR0[1]);
 477   fTPCinR0[4] = TMath::ATan2(fTPCinR0[1],fTPCinR0[0]);
 478   //
 479   fTPCinP0[0] = ref->Px();
 480   fTPCinP0[1] = ref->Py();
 481   fTPCinP0[2] = ref->Pz();
 482   fTPCinP0[3] = ref->Pt();
 483   fTPCinP0[4] = ref->P();
 484   fDeltaP     = (ref->P()-info->fParticle.P())/info->fParticle.P();
 485   //
 486   //
 487   if (fTPCinP0[3]>0.0000001){
 488     //
 489     fTPCAngle0[0] = TMath::ATan2(fTPCinP0[1],fTPCinP0[0]);
 490     fTPCAngle0[1] = TMath::ATan(fTPCinP0[2]/fTPCinP0[3]);
 491   }
 492   //
 493   if (fStatus[1]>0 &&info->fNTPCRef>0&&TMath::Abs(fTPCinP0[3])>0.0001){
 494     //TPC
 495     fESDtrack->GetInnerXYZ(fTPCinR1);
 496     fTPCinR1[3] = TMath::Sqrt(fTPCinR1[0]*fTPCinR1[0]+fTPCinR1[1]*fTPCinR1[1]);
 497     fTPCinR1[4] = TMath::ATan2(fTPCinR1[1],fTPCinR1[0]);
 498     fESDtrack->GetInnerPxPyPz(fTPCinP1);
 499     fTPCinP1[3] = TMath::Sqrt(fTPCinP1[0]*fTPCinP1[0]+fTPCinP1[1]*fTPCinP1[1]);
 500     fTPCinP1[4] = TMath::Sqrt(fTPCinP1[3]*fTPCinP1[3]+fTPCinP1[2]*fTPCinP1[2]);
 501     //
 502     //
 503     if (fTPCinP1[3]>0.000000000000001){
 504       fTPCAngle1[0] = TMath::ATan2(fTPCinP1[1],fTPCinP1[0]);
 505       fTPCAngle1[1] = TMath::ATan(fTPCinP1[2]/fTPCinP1[3]);
 506     }
 507     Double_t cov[15], param[5],x, alpha;
 508     fESDtrack->GetInnerExternalCovariance(cov);
 509     fESDtrack->GetInnerExternalParameters(alpha, x,param);
 510     if (x<50) return ;
 511     //
 512     fTPCDelta[0] = (fTPCinR0[4]-fTPCinR1[4])*fTPCinR1[3];  //delta rfi
 513     fTPCPools[0] = fTPCDelta[0]/TMath::Sqrt(cov[0]);
 514     fTPCDelta[1] = (fTPCinR0[2]-fTPCinR1[2]);              //delta z
 515     fTPCPools[1] = fTPCDelta[1]/TMath::Sqrt(cov[2]);
 516     fTPCDelta[2] = (fTPCAngle0[0]-fTPCAngle1[0]);
 517     fTPCPools[2] = fTPCDelta[2]/TMath::Sqrt(cov[5]);
 518     fTPCDelta[3] = (TMath::Tan(fTPCAngle0[1])-TMath::Tan(fTPCAngle1[1]));
 519     fTPCPools[3] = fTPCDelta[3]/TMath::Sqrt(cov[9]);
 520     fTPCDelta[4] = (fTPCinP0[3]-fTPCinP1[3]);
 521     Double_t sign = (param[4]>0)? 1.:-1;
 522     fSign =sign;
 523     fTPCPools[4] = sign*(1./fTPCinP0[3]-1./fTPCinP1[3])/TMath::Sqrt(TMath::Abs(cov[14]));
 524   }
 525 }
 526
 527
 528
 529
 530