HLT/comp/AliL3ModelTrack.cxx

   1 //$Id$
   2
   3 // Author: Anders Vestbo <mailto:vestbo$fi.uib.no>
   4 //*-- Copyright &copy ASV
   5
   6 #include "AliL3StandardIncludes.h"
   7
   8 #include "AliL3Logging.h"
   9 #include "AliL3ModelTrack.h"
  10 #include "AliL3Transform.h"
  11
  12 #if GCCVERSION == 3
  13 using namespace std;
  14 #endif
  15
  16 //_____________________________________________________________
  17 // AliL3ModelTrack
  18 //
  19 //
  20
  21 ClassImp(AliL3ModelTrack)
  22
  23 AliL3ModelTrack::AliL3ModelTrack()
  24 {
  25   fNClusters = 0;
  26   fClusters = 0;
  27   fOverlap = 0;
  28   fPad=0;
  29   fTime=0;
  30   fClusterCharge=0;
  31   fTrackModel=0;
  32   fLabel=0;
  33 }
  34
  35
  36 AliL3ModelTrack::~AliL3ModelTrack()
  37 {
  38   if(fClusters)
  39     delete [] fClusters;
  40   if(fPad)
  41     delete [] fPad;
  42   if(fTime)
  43     delete [] fTime;
  44   if(fOverlap)
  45     delete [] fOverlap;
  46   if(fTrackModel)
  47     delete fTrackModel;
  48 }
  49
  50 void AliL3ModelTrack::Init(Int_t slice,Int_t patch)
  51 {
  52   fNClusters = 0;
  53   fSlice=slice;
  54   fPatch=patch;
  55   Int_t nrows = AliL3Transform::GetNRows(fPatch);
  56   fClusters = new AliL3ClusterModel[nrows];
  57   fPad = new Float_t[nrows];
  58   fTime = new Float_t[nrows];
  59   fTrackModel = new AliL3TrackModel;
  60   fOverlap = new Int_t[nrows];
  61
  62   memset(fClusters,0,nrows*sizeof(AliL3ClusterModel));
  63   memset(fPad,0,nrows*sizeof(Float_t));
  64   memset(fTime,0,nrows*sizeof(Float_t));
  65   memset(fTrackModel,0,sizeof(AliL3TrackModel));
  66   for(Int_t i=0; i<nrows; i++)
  67     fOverlap[i]=-1;
  68 #ifdef do_mc
  69   for(Int_t i=0; i<nrows; i++)
  70     fClusters[i].fTrackID[0]=fClusters[i].fTrackID[1]=fClusters[i].fTrackID[2]=-2;
  71 #endif
  72   fClusterCharge = 100;
  73
  74   // 100 micrometers:
  75   fXYResidualQ = 0.01/AliL3Transform::GetPadPitchWidth(patch);
  76   fZResidualQ = 0.01/AliL3Transform::GetPadPitchWidth(patch);
  77
  78   fXYWidthQ = 0.005/AliL3Transform::GetPadPitchWidth(patch);
  79   fZWidthQ = 0.005/AliL3Transform::GetPadPitchWidth(patch);
  80 }
  81
  82
  83 void AliL3ModelTrack::SetCluster(Int_t row,Float_t fpad,Float_t ftime,Float_t charge,Float_t sigmaY2,Float_t sigmaZ2,Int_t npads)
  84 {
  85   Int_t index = row - AliL3Transform::GetFirstRow(fPatch);
  86   if(index != fNClusters)
  87     cout<<"AliL3ModelTrack::SetCluster() : Mismatch ; index: "<<index<<" nclusters "<<fNClusters<<endl;
  88
  89   if(index < 0 || index > AliL3Transform::GetNRows(fPatch))
  90     {
  91       cerr<<"AliL3ModelTrack::SetCluster() : Wrong index: "<<index<<" row "<<row<<endl;
  92       return;
  93     }
  94   AliL3ClusterModel *cl = GetClusterModel(row);
  95   if(!charge)
  96     cl->fPresent = kFALSE;
  97   else
  98     {
  99       cl->fPresent = kTRUE;
 100       cl->fDTime = (ftime - GetTimeHit(row))/fXYResidualQ;
 101       cl->fDPad = (fpad - GetPadHit(row))/fZResidualQ;
 102       cl->fDCharge = charge - fClusterCharge;
 103       cl->fDSigmaY2 = (sigmaY2 - GetParSigmaY2(row))/fXYWidthQ;
 104       cl->fDSigmaZ2 = (sigmaZ2 - GetParSigmaZ2(row))/fZWidthQ;
 105       cl->fNPads = npads;
 106     }
 107
 108   fNClusters++;
 109 }
 110
 111 Int_t AliL3ModelTrack::CheckClustersQuality(UInt_t npads)
 112 {
 113
 114   //Check the quality of clusters,- remove clusters with less than
 115   //npads.
 116   //Returns the number of good clusters left.
 117
 118   Int_t count=0;
 119
 120   for(Int_t i=AliL3Transform::GetFirstRow(fPatch); i<=AliL3Transform::GetLastRow(fPatch); i++)
 121     {
 122       AliL3ClusterModel *cl = GetClusterModel(i);
 123       if(cl->fNPads < npads)
 124         cl->fPresent = kFALSE;
 125       if(cl->fPresent)
 126         count++;
 127     }
 128
 129   return count;
 130 }
 131
 132 void AliL3ModelTrack::FillModel()
 133 {
 134   //Fill the track structure
 135
 136   if(!fTrackModel)
 137     {
 138       cerr<<"AliL3ModelTrack::FillModel() : No trackmodel "<<endl;
 139       return;
 140     }
 141   fTrackModel->fKappa = GetKappa();
 142   fTrackModel->fFirstPointX = GetFirstPointX();
 143   fTrackModel->fFirstPointY = GetFirstPointY();
 144   fTrackModel->fFirstPointZ = GetFirstPointZ();
 145   fTrackModel->fTgl = GetTgl();
 146   fTrackModel->fPsi = GetPsi();
 147   fTrackModel->fLength = (Short_t)GetLength();
 148   fTrackModel->fClusterCharge = fClusterCharge;
 149   fTrackModel->fNClusters = fNClusters;
 150
 151 }
 152
 153 void AliL3ModelTrack::FillTrack()
 154 {
 155   //Fill the track parameters from the structure.
 156
 157   if(!fTrackModel)
 158     {
 159       cerr<<"AliL3ModelTrack::FillTrack() : No data!!"<<endl;
 160       return;
 161     }
 162   SetKappa(fTrackModel->fKappa);
 163   SetCharge((-1*(Int_t)copysign(1.,GetKappa())));
 164   SetFirstPoint(fTrackModel->fFirstPointX,fTrackModel->fFirstPointY,fTrackModel->fFirstPointZ);
 165   SetTgl(fTrackModel->fTgl);
 166   SetPsi(fTrackModel->fPsi);
 167   SetLength(fTrackModel->fLength);
 168   fClusterCharge=fTrackModel->fClusterCharge;
 169   fNClusters = fTrackModel->fNClusters;
 170   SetPt((BFACT*AliL3Transform::GetBField())/fabs(GetKappa()));
 171
 172   CalculateHelix();
 173
 174   Float_t hit[3];
 175   Int_t sector,row;
 176   for(Int_t i=AliL3Transform::GetFirstRow(fPatch); i<=AliL3Transform::GetLastRow(fPatch); i++)
 177     {
 178       AliL3ClusterModel *cl = GetClusterModel(i);
 179       if(!cl) continue;
 180       GetCrossingPoint(i,hit);
 181       AliL3Transform::Slice2Sector(fSlice,i,sector,row);
 182       AliL3Transform::Local2Raw(hit,sector,row);
 183       SetPadHit(i,hit[1]);
 184       SetTimeHit(i,hit[2]);
 185     }
 186 }
 187
 188
 189 void AliL3ModelTrack::SetTrackID(Int_t row,Int_t *trackID)
 190 {
 191 #ifdef do_mc
 192   AliL3ClusterModel *cluster = GetClusterModel(row);
 193   cluster->fTrackID[0] = trackID[0];
 194   cluster->fTrackID[1] = trackID[1];
 195   cluster->fTrackID[2] = trackID[2];
 196   return;
 197 #endif
 198   cerr<<"AliL3ModelTrack::SetTrackID : Compile with do_mc flag"<<endl;
 199 }
 200
 201
 202 void AliL3ModelTrack::SetPadHit(Int_t row,Float_t pad)
 203 {
 204   Int_t index = row-AliL3Transform::GetFirstRow(fPatch);
 205   if(index < 0 || index > AliL3Transform::GetNRows(fPatch))
 206     {
 207       cerr<<"AliL3ModelTrack::SetPadHit() : Wrong index: "<<index<<endl;
 208       return;
 209     }
 210   fPad[index]=pad;
 211
 212 }
 213
 214 void AliL3ModelTrack::SetTimeHit(Int_t row,Float_t time)
 215 {
 216   Int_t index = row-AliL3Transform::GetFirstRow(fPatch);
 217   if(index < 0 || index > AliL3Transform::GetNRows(fPatch))
 218     {
 219       cerr<<"AliL3ModelTrack::SetTimeHit() : Wrong index: "<<index<<endl;
 220       return;
 221     }
 222   fTime[index]=time;
 223 }
 224
 225 void AliL3ModelTrack::SetOverlap(Int_t row,Int_t id)
 226 {
 227   Int_t index = row-AliL3Transform::GetFirstRow(fPatch);
 228   if(index < 0 || index > AliL3Transform::GetNRows(fPatch))
 229     {
 230       cerr<<"AliL3ModelTrack::SetOverlap() : Wrong index: "<<index<<endl;
 231       return;
 232     }
 233   fOverlap[index]=id;
 234 }
 235
 236
 237 Int_t AliL3ModelTrack::GetTrackID(Int_t row,Int_t index)
 238 {
 239
 240 #ifdef do_mc
 241   AliL3ClusterModel *cl = GetClusterModel(row);
 242   return cl->fTrackID[index];
 243 #endif
 244   cerr<<"AliL3ModelTrack::GetTrackID : Compile with do_mc flag"<<endl;
 245 }
 246
 247
 248 Int_t AliL3ModelTrack::GetNPads(Int_t row)
 249 {
 250   AliL3ClusterModel *cl = GetClusterModel(row);
 251   return cl->fNPads;
 252 }
 253
 254 Bool_t AliL3ModelTrack::GetPad(Int_t row,Float_t &pad)
 255 {
 256   //(ftime - GetTimeHit(fNClusters))/fXYResidualQ;
 257   //(fpad - GetPadHit(fNClusters))/fZResidualQ;
 258
 259   AliL3ClusterModel *cl = GetClusterModel(row);
 260   pad = cl->fDPad*fXYResidualQ + GetPadHit(row);
 261
 262   return (Bool_t)cl->fPresent;
 263 }
 264
 265 Bool_t AliL3ModelTrack::GetTime(Int_t row,Float_t &time)
 266 {
 267   AliL3ClusterModel *cl = GetClusterModel(row);
 268   time = cl->fDTime*fZResidualQ + GetTimeHit(row);
 269
 270   return (Bool_t)cl->fPresent;
 271 }
 272
 273 Bool_t AliL3ModelTrack::GetClusterCharge(Int_t row,Int_t &charge)
 274 {
 275   AliL3ClusterModel *cl = GetClusterModel(row);
 276   charge = (Int_t)cl->fDCharge + fClusterCharge;
 277
 278   return (Bool_t)cl->fPresent;
 279 }
 280
 281 Bool_t AliL3ModelTrack::GetXYWidth(Int_t row,Float_t &width)
 282 {
 283   AliL3ClusterModel *cl = GetClusterModel(row);
 284   width = cl->fDSigmaY2*fXYWidthQ + GetParSigmaY2(row);
 285
 286   return (Bool_t)cl->fPresent;
 287 }
 288
 289 Bool_t AliL3ModelTrack::GetZWidth(Int_t row,Float_t &width)
 290 {
 291   AliL3ClusterModel *cl = GetClusterModel(row);
 292   width = cl->fDSigmaZ2*fZWidthQ + GetParSigmaZ2(row);
 293
 294   return (Bool_t)cl->fPresent;
 295 }
 296
 297 Bool_t AliL3ModelTrack::GetPadResidual(Int_t row,Float_t &res)
 298 {
 299   AliL3ClusterModel *cl = GetClusterModel(row);
 300   res = cl->fDPad;
 301   return cl->fPresent;
 302 }
 303
 304 Bool_t AliL3ModelTrack::GetTimeResidual(Int_t row,Float_t &res)
 305 {
 306   AliL3ClusterModel *cl = GetClusterModel(row);
 307   res = cl->fDTime;
 308   return cl->fPresent;
 309 }
 310
 311 Bool_t AliL3ModelTrack::GetXYWidthResidual(Int_t row,Float_t &res)
 312 {
 313   AliL3ClusterModel *cl = GetClusterModel(row);
 314   res = cl->fDSigmaY2;
 315   return cl->fPresent;
 316 }
 317
 318 Bool_t AliL3ModelTrack::GetZWidthResidual(Int_t row,Float_t &res)
 319 {
 320   AliL3ClusterModel *cl = GetClusterModel(row);
 321   res = cl->fDSigmaZ2;
 322   return cl->fPresent;
 323 }
 324
 325 Float_t AliL3ModelTrack::GetPadHit(Int_t row)
 326 {
 327   Int_t index = row-AliL3Transform::GetFirstRow(fPatch);
 328   if(index < 0 || index > AliL3Transform::GetNRows(fPatch))
 329     {
 330       cerr<<"AliL3ModelTrack::GetPadHit() : Wrong index: "<<index<<" row "<<row<<endl;
 331       return 0;
 332     }
 333   return fPad[index];
 334 }
 335
 336 Float_t AliL3ModelTrack::GetTimeHit(Int_t row)
 337 {
 338   Int_t index = row-AliL3Transform::GetFirstRow(fPatch);
 339   if(index < 0 || index > AliL3Transform::GetNRows(fPatch))
 340     {
 341       cerr<<"AliL3ModelTrack::GetTimeHit() : Wrong index: "<<index<<" row "<<row<<endl;
 342       return 0;
 343     }
 344   return fTime[index];
 345 }
 346
 347 Int_t AliL3ModelTrack::GetOverlap(Int_t row)
 348 {
 349   Int_t index = row-AliL3Transform::GetFirstRow(fPatch);
 350   if(index < 0 || index > AliL3Transform::GetNRows(fPatch))
 351     {
 352       cerr<<"AliL3ModelTrack::GetOverlap() : Wrong index: "<<index<<endl;
 353       return 0;
 354     }
 355   return fOverlap[index];
 356 }
 357
 358 AliL3ClusterModel *AliL3ModelTrack::GetClusterModel(Int_t row)
 359 {
 360   if(!fClusters) return 0;
 361   Int_t index = row-AliL3Transform::GetFirstRow(fPatch);
 362   if(index < 0 || index > AliL3Transform::GetNRows(fPatch))
 363     {
 364       cerr<<"AliL3ModelTrack::GetClusterModel() : Wrong index: "<<index<<endl;
 365       return 0;
 366     }
 367   return &fClusters[index];
 368 }
 369
 370 void AliL3ModelTrack::Print()
 371 {
 372   //Print info
 373
 374   cout<<"----Slice "<<fSlice<<" Patch "<<fPatch<<"----"<<endl;
 375   cout<<"First point "<<GetFirstPointX()<<" "<<GetFirstPointY()<<" "<<GetFirstPointZ()<<endl;
 376   cout<<"Last point "<<GetLastPointX()<<" "<<GetLastPointY()<<" "<<GetLastPointZ()<<endl;
 377   cout<<"Pt "<<GetPt()<<" kappa "<<GetKappa()<<" tgl "<<GetTgl()<<" psi "<<GetPsi()<<" charge "<<GetCharge()<<endl;
 378   cout<<"Center "<<GetCenterX()<<" "<<GetCenterY()<<endl<<endl;
 379   cout<<"NHits "<<GetNClusters()<<endl;
 380   cout<<"Clusters:"<<endl;
 381
 382   for(Int_t i=AliL3Transform::GetFirstRow(fPatch); i<=AliL3Transform::GetLastRow(fPatch); i++)
 383     {
 384       AliL3ClusterModel *cl = GetClusterModel(i);
 385
 386       if(!cl->fPresent)
 387         cout<<i<<" Empty"<<" Padcrossing "<<GetPadHit(i)<<" Timecrossing "<<GetTimeHit(i)<<" ";
 388       else
 389         {
 390           cout<<i<<" Dpad "<<cl->fDPad<<" Dtime "<<cl->fDTime<<" Dcharge "<<cl->fDCharge;
 391           cout<<" DsigmaY2 "<<cl->fDSigmaY2<<" DsigmaZ2 "<<cl->fDSigmaZ2;
 392           cout<<" Padcrossing "<<GetPadHit(i)<<" Timecrossing "<<GetTimeHit(i)<<" ";
 393           cout<<"Number of pads "<<GetNPads(i)<<endl;
 394         }
 395       cout<<"Overlapping index "<<GetOverlap(i)<<endl;
 396     }
 397 }
 398
 399 Double_t AliL3ModelTrack::GetParSigmaY2(Int_t row)
 400 {
 401   //Calculate the expected cluster width, based on the track parameters and drift distance.
 402
 403   Float_t pad,time;
 404   if(!GetTime(row,time) || !GetPad(row,pad))
 405     return -1;
 406
 407   Float_t xyz[3];
 408   Int_t sector,padrow;
 409   AliL3Transform::Slice2Sector(fSlice,row,sector,padrow);
 410   AliL3Transform::Raw2Local(xyz,sector,padrow,pad,time);
 411
 412   //Calculate the drift length:
 413   Double_t drift;
 414   if(xyz[2] > 0)
 415     drift = AliL3Transform::GetZLength() - xyz[2];
 416   else
 417     drift = AliL3Transform::GetZLength() + xyz[2];
 418
 419   Double_t prf = AliL3Transform::GetPRFSigma(fPatch);
 420   Double_t diffT = AliL3Transform::GetDiffT();
 421   Double_t padlength = AliL3Transform::GetPadLength(fPatch);
 422   Double_t anode = AliL3Transform::GetAnodeWireSpacing();
 423   Double_t beta = GetCrossingAngle(row);
 424
 425   Double_t sigmaY2 = prf*prf + diffT*diffT*drift + padlength*padlength*tan(beta)*tan(beta)/12 + anode*anode*pow( tan(beta)-0.15, 2)/12;
 426
 427   //Convert back to raw coordinates.
 428   sigmaY2 = sigmaY2/pow(AliL3Transform::GetPadPitchWidth(fPatch),2);
 429   return sigmaY2;
 430 }
 431
 432 Double_t AliL3ModelTrack::GetParSigmaZ2(Int_t row)
 433 {
 434   //Calculate the expected cluster width, based on the track parameters and drift distance.
 435
 436   Float_t pad,time;
 437   if(!GetTime(row,time) || !GetPad(row,pad))
 438     return -1;
 439
 440   Float_t xyz[3];
 441   Int_t sector,padrow;
 442   AliL3Transform::Slice2Sector(fSlice,row,sector,padrow);
 443   AliL3Transform::Raw2Local(xyz,sector,padrow,pad,time);
 444
 445   //Calculate the drift length:
 446   Double_t drift;
 447   if(xyz[2] > 0)
 448     drift = AliL3Transform::GetZLength() - xyz[2];
 449   else
 450     drift = AliL3Transform::GetZLength() + xyz[2];
 451
 452   Double_t sigma0 = AliL3Transform::GetTimeSigma();
 453   Double_t diffL = AliL3Transform::GetDiffL();
 454   Double_t padlength = AliL3Transform::GetPadLength(fPatch);
 455   Double_t tanl = GetTgl();
 456
 457   Double_t sigmaZ2 = sigma0*sigma0 + diffL*diffL*drift + padlength*padlength * tanl*tanl/12;
 458
 459   //Convert back to raw coodinates:
 460   sigmaZ2 = sigmaZ2/pow(AliL3Transform::GetZWidth(),2);
 461   return sigmaZ2;
 462
 463 }
 464
 465 void AliL3ModelTrack::AssignTrackID(Float_t wrong)
 466 {
 467   //Assign a track ID to the track, corresponding to the MC TParticle ID.
 468   //Can only be done if you compiled with do_mc flag, of course.
 469   //The function loops over the assigned clusters, and finds the label (ID)
 470   //of each clusters, and assigns the ID with the most hits to the track.
 471   //If there are more than wrong% clusters of a different ID, the track is
 472   //considered to be fake, and label will be assigned as negative.
 473
 474 #ifdef do_mc
 475   Int_t *lb = new Int_t[GetNClusters()];
 476   Int_t *mx = new Int_t[GetNClusters()];
 477
 478   Int_t i,j;
 479   for(Int_t i=0; i<GetNClusters(); i++)
 480     lb[i]=mx[i]=0;
 481
 482   Int_t lab=123456789;
 483
 484   for(i=0; i<GetNClusters(); i++)
 485     {
 486       lab = abs(GetTrackID(i,0));
 487       for (j=0; j<GetNClusters(); j++)
 488         if (lb[j]==lab || mx[j]==0) break;
 489       lb[j]=lab;
 490       (mx[j])++;
 491     }
 492
 493   Int_t max=0;
 494   for (i=0; i<GetNClusters(); i++)
 495     {
 496       if(mx[i] > max)
 497         {
 498           max=mx[i];
 499           lab=lb[i];
 500         }
 501     }
 502
 503   for (i=0; i<GetNClusters(); i++)
 504     {
 505       if(abs(GetTrackID(i,1)) == lab ||
 506          abs(GetTrackID(i,2)) == lab)
 507         max++;
 508     }
 509
 510   if ((1.- Float_t(max)/GetNClusters()) > wrong) lab=-lab;
 511
 512   SetLabel(lab);
 513
 514   delete[] lb;
 515   delete[] mx;
 516   return;
 517 #endif
 518   cerr<<"AliL3ModelTrack::AssignTrackID : Compile with do_mc flag"<<endl;
 519 }