TOF/AliTOFtrackerV2.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 // AliTOFtrackerV2 Class                                              //
  19 // Task: Perform association of the ESD tracks to TOF Clusters        //
  20 // and Update ESD track with associated TOF Cluster parameters        //
  21 //                                                                    //
  22 // -- Authors : S. Arcelli, C. Zampolli (Bologna University and INFN) //
  23 // -- Contacts: Annalisa.De.Caro@cern.ch                              //
  24 // --         : Chiara.Zampolli@bo.infn.it                            //
  25 // --         : Silvia.Arcelli@bo.infn.it                             //
  26 //                                                                    //
  27 //--------------------------------------------------------------------//
  28
  29 #include <Rtypes.h>
  30 #include <TROOT.h>
  31
  32 #include <TClonesArray.h>
  33 #include <TObjArray.h>
  34 #include <TGeoManager.h>
  35 #include <TTree.h>
  36
  37 #include "AliGeomManager.h"
  38 #include "AliESDtrack.h"
  39 #include "AliESDEvent.h"
  40 #include "AliESDpid.h"
  41 #include "AliESDTOFCluster.h"
  42 #include "AliLog.h"
  43 #include "AliTrackPointArray.h"
  44 #include "AliCDBManager.h"
  45
  46 #include "AliTOFRecoParam.h"
  47 #include "AliTOFReconstructor.h"
  48 #include "AliTOFGeometry.h"
  49 #include "AliTOFtrackerV2.h"
  50 #include "AliTOFtrack.h"
  51
  52 #include "AliESDTOFHit.h"
  53
  54 extern TGeoManager *gGeoManager;
  55
  56 ClassImp(AliTOFtrackerV2)
  57
  58 //_____________________________________________________________________________
  59 AliTOFtrackerV2::AliTOFtrackerV2():
  60   fkRecoParam(0x0),
  61   fN(0),
  62   fNseeds(0),
  63   fNseedsTOF(0),
  64   fnunmatch(0),
  65   fnmatch(0),
  66   fSeeds(new TObjArray(100)),
  67   fClustersESD(new TClonesArray("AliESDTOFCluster")),
  68   fHitsESD(new TClonesArray("AliESDTOFHit")),
  69   fEvent(0),
  70   fNsteps(0)
  71 {
  72   //AliTOFtrackerV2 main Ctor
  73   for (Int_t ii=0; ii<kMaxCluster; ii++){
  74     fClusters[ii]=0x0;
  75     fWrittenInPos[ii] = -1;
  76   }
  77
  78   for(Int_t isp=0;isp < AliPID::kSPECIESC;isp++)
  79     fTimesAr[isp] = NULL;
  80
  81   for (Int_t ii=0; ii<4; ii++)
  82     fTrackPos[ii] = NULL;
  83 }
  84 //_____________________________________________________________________________
  85 AliTOFtrackerV2::~AliTOFtrackerV2() {
  86   //
  87   // Dtor
  88   //
  89
  90   if(!(AliCDBManager::Instance()->GetCacheFlag())){
  91     delete fkRecoParam;
  92   }
  93   if (fSeeds){
  94     fSeeds->Delete();
  95     delete fSeeds;
  96     fSeeds=0x0;
  97   }
  98
  99   if (fClustersESD){
 100     fClustersESD->Delete();
 101     delete fClustersESD;
 102     fClustersESD=0x0;
 103   }
 104   if (fHitsESD){
 105     fHitsESD->Delete();
 106     delete fHitsESD;
 107     fHitsESD=0x0;
 108   }
 109
 110   for(Int_t isp=0;isp < AliPID::kSPECIESC;isp++){
 111     if(fTimesAr[isp]) delete[] fTimesAr[isp];
 112     fTimesAr[isp] = NULL;
 113   }
 114
 115
 116   for (Int_t ii=0; ii<4; ii++)
 117     if(fTrackPos[ii])
 118       delete [] fTrackPos[ii];
 119 }
 120 //_____________________________________________________________________________
 121 void AliTOFtrackerV2::GetPidSettings(AliESDpid *esdPID) {
 122   //
 123   // Sets TOF resolution from RecoParams
 124   //
 125   if (fkRecoParam)
 126     esdPID->GetTOFResponse().SetTimeResolution(fkRecoParam->GetTimeResolution());
 127   else
 128     AliWarning("fkRecoParam not yet set; cannot set PID settings");
 129 }
 130 //_____________________________________________________________________________
 131 Int_t AliTOFtrackerV2::PropagateBack(AliESDEvent * const event) {
 132   //
 133   // Gets seeds from ESD event and Match with TOF Clusters
 134   //
 135
 136   //Update the matched ESD tracks
 137   // needed in case of call of TOF info before of the selection of matching and in case of no clusters available at all
 138
 139   fEvent = event;
 140
 141   if (fN==0) {
 142     AliInfo("No TOF recPoints to be matched with reconstructed tracks");
 143     return 0;
 144   }
 145
 146   // initialize RecoParam for current event
 147   AliDebug(1,"Initializing params for TOF");
 148
 149   fkRecoParam = AliTOFReconstructor::GetRecoParam();  // instantiate reco param from STEER...
 150
 151   if (fkRecoParam == 0x0) {
 152     AliFatal("No Reco Param found for TOF!!!");
 153   }
 154
 155   //Initialise some counters
 156
 157   fNseeds=0;
 158   fNseedsTOF=0;
 159
 160   Int_t ntrk=event->GetNumberOfTracks();
 161   fNseeds = ntrk;
 162
 163   //Load ESD tracks into a local Array of ESD Seeds
 164   for (Int_t i=0; i<fNseeds; i++){
 165     fSeeds->AddLast(event->GetTrack(i));
 166     event->GetTrack(i)->SetESDEvent(event);
 167   }
 168   //Prepare ESD tracks candidates for TOF Matching
 169   CollectESD();
 170
 171   if (fNseeds==0 || fNseedsTOF==0) {
 172     AliInfo("No seeds to try TOF match");
 173     fSeeds->Clear();
 174     fClustersESD->Clear();
 175     fHitsESD->Clear();
 176     return 0;
 177   }
 178
 179   // clusterize before of matching
 180   Clusterize(); // fN might change
 181
 182   //Second Step with Looser Matching Criterion
 183   MatchTracks();
 184
 185   // switch array from ALL to filter for ESD (moving from all fClusterESD to filtered AliESDEvent->GetESDTOFClusters())
 186   TClonesArray* esdTOFHitArr = event->GetESDTOFHits();
 187   TClonesArray* esdTOFClArr = event->GetESDTOFClusters();
 188
 189   AliInfo(Form("TOF before the matching: hits = %i and clusters = %i",esdTOFHitArr->GetEntriesFast(),fClustersESD->GetEntriesFast()));
 190
 191   while(esdTOFHitArr->GetEntriesFast()){ // remove all hits
 192     delete esdTOFHitArr->RemoveAt(esdTOFHitArr->GetEntriesFast()-1);
 193   }
 194   for(Int_t i=0;i < fHitsESD->GetEntriesFast();i++){
 195     AliESDTOFHit *hitToStored = (AliESDTOFHit *) fHitsESD->At(i);
 196     AliESDTOFHit *hitNew = new ( (*esdTOFHitArr)[esdTOFHitArr->GetEntriesFast()] ) AliESDTOFHit(*hitToStored);
 197     hitNew->SetESDTOFClusterIndex(hitToStored->GetESDTOFClusterIndex());
 198   }
 199
 200   AliInfo(Form("TOF after the matching: hits = %i and clusters = %i",esdTOFHitArr->GetEntriesFast(),esdTOFClArr->GetEntriesFast()));
 201
 202   // Sort tof cluster
 203   for (Int_t i=0; i<fNseeds; i++) {
 204     AliESDtrack *t =(AliESDtrack*)fSeeds->At(i);
 205     if ((t->GetStatus()&AliESDtrack::kTOFout)==0)continue;
 206     t->SortTOFcluster();
 207   }
 208
 209
 210   fSeeds->Clear();
 211   fClustersESD->Clear();
 212   fHitsESD->Clear();
 213   return 0;
 214
 215 }
 216 //_________________________________________________________________________
 217 void AliTOFtrackerV2::CollectESD() {
 218    //prepare the set of ESD tracks to be matched to clusters in TOF
 219
 220   Int_t seedsTOF1=0;
 221   Int_t seedsTOF3=0;
 222   Int_t seedsTOF2=0;
 223
 224   AliTOFtrack track;
 225
 226   for (Int_t i=0; i<fNseeds; i++) {
 227
 228     AliESDtrack *t =(AliESDtrack*)fSeeds->At(i);
 229     if ((t->GetStatus()&AliESDtrack::kTPCout)==0)continue;
 230
 231     track = *t; // New
 232     Float_t x = (Float_t)track.GetX(); //New
 233
 234     // TRD 'good' tracks
 235     if ( ( (t->GetStatus()&AliESDtrack::kTRDout)!=0 ) ) {
 236
 237       AliDebug(1,Form(" Before propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track.GetIntegratedLength()));
 238
 239       // TRD 'good' tracks, already propagated at 371 cm
 240       if ( x >= AliTOFGeometry::Rmin() ) {
 241
 242         if ( track.PropagateToInnerTOF() ) {
 243
 244           AliDebug(1,Form(" TRD propagated track till rho = %fcm."
 245                           " And then the track has been propagated till rho = %fcm.",
 246                           x, (Float_t)track.GetX()));
 247
 248           track.SetSeedIndex(i);
 249           t->UpdateTrackParams(&track,AliESDtrack::kTOFin);
 250           fNseedsTOF++;
 251           seedsTOF1++;
 252
 253           AliDebug(1,Form(" After propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track.GetIntegratedLength()));
 254         }
 255       }
 256       else { // TRD 'good' tracks, propagated rho<371cm
 257
 258         if  ( track.PropagateToInnerTOF() ) {
 259
 260           AliDebug(1,Form(" TRD propagated track till rho = %fcm."
 261                           " And then the track has been propagated till rho = %fcm.",
 262                           x, (Float_t)track.GetX()));
 263
 264           track.SetSeedIndex(i);
 265           t->UpdateTrackParams(&track,AliESDtrack::kTOFin);
 266           fNseedsTOF++;
 267           seedsTOF3++;
 268
 269           AliDebug(1,Form(" After propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track.GetIntegratedLength()));
 270         }
 271       }
 272     }
 273
 274     else { // Propagate the rest of TPCbp
 275
 276       AliDebug(1,Form(" Before propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track.GetIntegratedLength()));
 277
 278       if ( track.PropagateToInnerTOF() ) {
 279
 280         AliDebug(1,Form(" TPC propagated track till rho = %fcm."
 281                         " And then the track has been propagated till rho = %fcm.",
 282                         x, (Float_t)track.GetX()));
 283
 284         track.SetSeedIndex(i);
 285         t->UpdateTrackParams(&track,AliESDtrack::kTOFin);
 286         fNseedsTOF++;
 287         seedsTOF2++;
 288
 289         AliDebug(1,Form(" After propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track.GetIntegratedLength()));
 290       }
 291     }
 292   }
 293
 294   AliInfo(Form("Number of TOF seeds = %d (kTRDout371 = %d, kTRDoutLess371 = %d, !kTRDout = %d)",fNseedsTOF,seedsTOF1,seedsTOF3,seedsTOF2));
 295
 296 }
 297
 298 //_________________________________________________________________________
 299 void AliTOFtrackerV2::MatchTracks() {
 300   //
 301   //Match ESD tracks to clusters in TOF
 302   //
 303
 304   // Parameters used/regulating the reconstruction
 305   static Float_t detDepth=18.;
 306   static Float_t padDepth=0.5;
 307
 308   const Float_t kSpeedOfLight= 2.99792458e-2; // speed of light [cm/ps]
 309
 310   Float_t dY=AliTOFGeometry::XPad();
 311   Float_t dZ=AliTOFGeometry::ZPad();
 312
 313   Float_t sensRadius = fkRecoParam->GetSensRadius();
 314   Float_t stepSize   = fkRecoParam->GetStepSize();
 315   Float_t scaleFact  = fkRecoParam->GetWindowScaleFact();
 316   Float_t dyMax=fkRecoParam->GetWindowSizeMaxY();
 317   Float_t dzMax=fkRecoParam->GetWindowSizeMaxZ();
 318   Float_t dCut=10.;//fkRecoParam->GetDistanceCut(); // This is to be loaded by OCDB. It should be 10cm always.
 319   Double_t maxChi2=fkRecoParam->GetMaxChi2TRD();
 320   Bool_t timeWalkCorr = fkRecoParam->GetTimeWalkCorr();
 321   AliDebug(1,"++++++++++++++TOF Reconstruction Parameters:++++++++++++++");
 322   AliDebug(1,Form("TOF sens radius: %f",sensRadius));
 323   AliDebug(1,Form("TOF step size: %f",stepSize));
 324   AliDebug(1,Form("TOF Window scale factor: %f",scaleFact));
 325   AliDebug(1,Form("TOF Window max dy: %f",dyMax));
 326   AliDebug(1,Form("TOF Window max dz: %f",dzMax));
 327   AliDebug(1,Form("TOF distance Cut: %f",dCut));
 328   AliDebug(1,Form("TOF Max Chi2: %f",maxChi2));
 329   AliDebug(1,Form("Time Walk Correction? : %d",timeWalkCorr));
 330
 331   //Match ESD tracks to clusters in TOF
 332   TClonesArray* TOFClArr = fClustersESD; // use temporary array
 333   TClonesArray* esdTOFClArr = fEvent->GetESDTOFClusters();
 334   TClonesArray* esdTOFHitArr = fEvent->GetESDTOFHits();
 335
 336   if(Int_t(detDepth/stepSize) > fNsteps){ // create array for each step
 337     // Get the number of propagation steps
 338     fNsteps =(Int_t)(detDepth/stepSize);
 339
 340     for(Int_t isp=0;isp < AliPID::kSPECIESC;isp++){
 341       if(fTimesAr[isp]) delete[] fTimesAr[isp];
 342     }
 343
 344     for(Int_t isp=0;isp < AliPID::kSPECIESC;isp++){
 345       fTimesAr[isp] = new Double_t[fNsteps];
 346     }
 347
 348     for (Int_t ii=0; ii<4; ii++)
 349       if(fTrackPos[ii])
 350         delete [] fTrackPos[ii];
 351
 352     for (Int_t ii=0; ii<4; ii++) fTrackPos[ii] = new Float_t[fNsteps];
 353   }
 354
 355
 356
 357   AliDebug(1,Form(" Number of steps to be done %d",fNsteps));
 358
 359   AliDebug(1,"++++++++++++++++++++++++++++++++++++++++++++++++++++++++++");
 360
 361   // Some init
 362   const Int_t kNclusterMax = 1000; // related to fN value
 363   TGeoHMatrix global[kNclusterMax];
 364   Int_t clind[kNclusterMax];
 365   Bool_t isClusterMatchable[kNclusterMax]; // true if track and cluster were already matched (set to false below upto nc < kNclusterMax)
 366
 367   AliTOFtrack trackTOFin;
 368
 369   //The matching loop
 370   for (Int_t iseed=0; iseed<fSeeds->GetEntriesFast(); iseed++) {
 371     AliESDtrack *t =(AliESDtrack*)fSeeds->At(iseed); // ciao replace with loop on ESD + kTOFin
 372     if( (t->GetStatus()&AliESDtrack::kTOFin) == 0 ) continue;
 373
 374     trackTOFin = *t;
 375
 376     for (Int_t ii=0; ii<4; ii++)
 377       for (Int_t jj=0; jj<fNsteps; jj++) fTrackPos[ii][jj]=0.;
 378
 379     for (Int_t ii=0; ii<kNclusterMax; ii++) clind[ii]=-1;
 380     for (Int_t ii=0; ii<kNclusterMax; ii++) global[ii] = 0x0;
 381     for (Int_t ii=0; ii<kNclusterMax; ii++) isClusterMatchable[ii] = kFALSE;
 382
 383     Double_t timesOr[AliPID::kSPECIESC]; t->GetIntegratedTimes(timesOr,AliPID::kSPECIESC); // in ps
 384
 385     // Determine a window around the track
 386     Double_t x,par[5];
 387     trackTOFin.GetExternalParameters(x,par);
 388     Double_t cov[15];
 389     trackTOFin.GetExternalCovariance(cov);
 390
 391     if (cov[0]<0. || cov[2]<0.) {
 392       AliWarning(Form("Very strange track (%d)! At least one of its covariance matrix diagonal elements is negative!",iseed));
 393       continue;
 394     }
 395
 396     Double_t dphi=
 397       scaleFact*
 398       ((5*TMath::Sqrt(TMath::Abs(cov[0])) + 0.5*dY + 2.5*TMath::Abs(par[2]))/sensRadius);
 399     Double_t dz=
 400        scaleFact*
 401        (5*TMath::Sqrt(TMath::Abs(cov[2])) + 0.5*dZ + 2.5*TMath::Abs(par[3]));
 402
 403     Double_t phi=TMath::ATan2(par[0],x) + trackTOFin.GetAlpha();
 404     if (phi<-TMath::Pi())phi+=2*TMath::Pi();
 405     if (phi>=TMath::Pi())phi-=2*TMath::Pi();
 406     Double_t z=par[1];
 407
 408     //upper limit on window's size.
 409     if (dz> dzMax) dz=dzMax;
 410     if (dphi*sensRadius> dyMax) dphi=dyMax/sensRadius;
 411
 412
 413     // find the clusters in the window of the track
 414     Int_t nc=0;
 415     for (Int_t k=FindClusterIndex(z-dz); k<TOFClArr->GetEntriesFast(); k++) {
 416
 417       if (nc>=kNclusterMax) {
 418         AliWarning("No more matchable clusters can be stored! Please, increase the corresponding vectors size.");
 419         break;
 420       }
 421
 422       AliESDTOFCluster *c=(AliESDTOFCluster *) TOFClArr->At(k);
 423       if (c->GetZ() > z+dz) break;
 424       if (!c->GetStatus()) {
 425         AliDebug(1,"Cluster in channel declared bad!");
 426         continue; // skip bad channels as declared in OCDB
 427       }
 428
 429
 430       Double_t dph=TMath::Abs(c->GetPhi()-phi);
 431       if (dph>TMath::Pi()) dph-=2.*TMath::Pi();
 432       if (TMath::Abs(dph)>dphi) continue;
 433
 434       Double_t yc=(c->GetPhi() - trackTOFin.GetAlpha())*c->GetR();
 435       Double_t p[2]={yc, c->GetZ()};
 436       Double_t cov2[3]= {dY*dY/12., 0., dZ*dZ/12.};
 437       if (trackTOFin.AliExternalTrackParam::GetPredictedChi2(p,cov2) > maxChi2)continue;
 438
 439       clind[nc] = k;
 440       Char_t path[200];
 441       Int_t ind[5]; AliTOFGeometry::GetVolumeIndices(c->GetTOFchannel(),ind);
 442       AliTOFGeometry::GetVolumePath(ind,path);
 443       gGeoManager->cd(path);
 444       global[nc] = *gGeoManager->GetCurrentMatrix();
 445       nc++;
 446     }
 447
 448
 449     if (nc == 0 ) {
 450       AliDebug(1,Form("No available clusters for the track number %d",iseed));
 451       fnunmatch++;
 452       continue;
 453     }
 454
 455     AliDebug(1,Form(" Number of available TOF clusters for the track number %d: %d",iseed,nc));
 456
 457     //start fine propagation
 458
 459     Int_t nStepsDone = 0;
 460     for( Int_t istep=0; istep<fNsteps; istep++){
 461
 462       // First of all, propagate the track...
 463       Float_t xs = AliTOFGeometry::RinTOF()+istep*stepSize;
 464       if (!(trackTOFin.PropagateTo(xs))) break;
 465
 466       //  ...and then, if necessary, rotate the track
 467       Double_t ymax = xs*TMath::Tan(0.5*AliTOFGeometry::GetAlpha());
 468       Double_t ysect = trackTOFin.GetY();
 469       if (ysect > ymax) {
 470         if (!(trackTOFin.Rotate(AliTOFGeometry::GetAlpha()))) break;
 471       } else if (ysect <-ymax) {
 472         if (!(trackTOFin.Rotate(-AliTOFGeometry::GetAlpha()))) break;
 473       }
 474
 475       Double_t mom = trackTOFin.P();
 476
 477       if(istep == 0){
 478         for(Int_t isp=0;isp<AliPID::kSPECIESC;isp++){
 479           Double_t mass=AliPID::ParticleMass(isp);
 480           Double_t momz = mom*AliPID::ParticleCharge(isp);
 481           fTimesAr[isp][nStepsDone] = stepSize/kSpeedOfLight*TMath::Sqrt(momz*momz+mass*mass)/momz;
 482         }
 483       }
 484       else{
 485         for(Int_t isp=0;isp<AliPID::kSPECIESC;isp++){
 486           Double_t mass=AliPID::ParticleMass(isp);
 487           Double_t momz = mom*AliPID::ParticleCharge(isp);
 488           fTimesAr[isp][nStepsDone] = fTimesAr[isp][nStepsDone-1] + (trackTOFin.GetIntegratedLength()-fTrackPos[3][nStepsDone-1])/kSpeedOfLight*TMath::Sqrt(momz*momz+mass*mass)/momz;
 489         }
 490       }
 491
 492       // store the running point (Globalrf) - fine propagation
 493
 494       Double_t r[3]; trackTOFin.GetXYZ(r);
 495       fTrackPos[0][nStepsDone]= (Float_t) r[0];
 496       fTrackPos[1][nStepsDone]= (Float_t) r[1];
 497       fTrackPos[2][nStepsDone]= (Float_t) r[2];
 498       fTrackPos[3][nStepsDone]= trackTOFin.GetIntegratedLength();
 499
 500       nStepsDone++;
 501       AliDebug(3,Form(" current step %d (%d) - nStepsDone=%d",istep,fNsteps,nStepsDone));
 502     }
 503
 504     if ( nStepsDone == 0 ) {
 505       AliDebug(1,Form(" No track points for track number %d",iseed));
 506       fnunmatch++;
 507       continue;
 508     }
 509
 510     AliDebug(3,Form(" Number of steps done for the track number %d: %d",iseed,nStepsDone));
 511
 512     if(nc){
 513       for (Int_t i=0; i<nc; i++) isClusterMatchable[i] = kFALSE;
 514     }
 515
 516     Int_t nfound = 0;
 517     Bool_t accept = kFALSE;
 518     for (Int_t istep=0; istep<nStepsDone; istep++) {
 519       Float_t ctrackPos[3];
 520       ctrackPos[0] = fTrackPos[0][istep];
 521       ctrackPos[1] = fTrackPos[1][istep];
 522       ctrackPos[2] = fTrackPos[2][istep];
 523
 524       //now see whether the track matches any of the TOF clusters
 525
 526       Float_t dist3d[3]={0.,0.,0.};
 527       accept = kFALSE;
 528
 529       for (Int_t i=0; i<nc; i++) {
 530
 531         AliTOFGeometry::IsInsideThePad((TGeoHMatrix*)(&global[i]),ctrackPos,dist3d);
 532
 533         // check multiple hit cases
 534         AliESDTOFCluster *cmatched=(AliESDTOFCluster *) TOFClArr->At(clind[i]);
 535
 536         if(cmatched->GetNTOFhits() > 1){ // correct residual for mean position of the clusters (w.r.t. the first pad/hit)
 537           Float_t zmain = cmatched->GetTOFchannel(0)/48;
 538           Float_t xmain = cmatched->GetTOFchannel(0)%48;
 539           for(Int_t ihit=1;ihit < cmatched->GetNTOFhits();ihit++){
 540             Float_t deltaz = (cmatched->GetTOFchannel(ihit)/48 - zmain) * 3.5;
 541             Float_t deltax = (cmatched->GetTOFchannel(ihit)%48 - xmain) * 2.5;
 542             dist3d[0] -= deltax / cmatched->GetNTOFhits();
 543             dist3d[2] -= deltaz / cmatched->GetNTOFhits();
 544           }
 545         }
 546
 547         // ***** NEW *****
 548         /* if track is inside this cluster set flags which will then
 549          * inhibit to add track points for the other clusters */
 550         Float_t yLoc = dist3d[1];
 551         Float_t rLoc = TMath::Sqrt(dist3d[0]*dist3d[0]+dist3d[2]*dist3d[2]);
 552         accept = (TMath::Abs(yLoc)<padDepth*0.5 && rLoc<dCut);
 553
 554         //***** NEW *****
 555         /* add point everytime that:
 556          * - the tracks is within dCut from the cluster
 557          */
 558         if (accept) {
 559
 560           Double_t timesCurrent[AliPID::kSPECIESC];
 561           AliDebug(3,Form(" Momentum for track %d -> %f", iseed,t->P()));
 562           for (Int_t j=0;j<AliPID::kSPECIESC;j++) {
 563             timesCurrent[j] = timesOr[j] + fTimesAr[j][istep];
 564           }
 565
 566
 567           if (TMath::Abs(dist3d[1])<stepSize && !isClusterMatchable[i]) {
 568             isClusterMatchable[i] = kTRUE;
 569
 570             Int_t currentpos = esdTOFClArr->GetEntriesFast(); // position of cluster in ESD
 571             if(fWrittenInPos[clind[i]] != -1){
 572               currentpos = fWrittenInPos[clind[i]];
 573               cmatched = (AliESDTOFCluster *) esdTOFClArr->At(currentpos); // update the new one in the ESDEvent
 574             }
 575             else{ // add as a new cluster in the ESD TClonesArray
 576               AliESDTOFCluster *clnew =  new( (*esdTOFClArr)[currentpos] ) AliESDTOFCluster(*cmatched);
 577               clnew->SetEvent(fEvent);
 578               clnew->SetESDID(currentpos);
 579
 580               // remap also TOF hit in the filtered array
 581               for(Int_t ii=0;ii < cmatched->GetNTOFhits();ii++){
 582                 Int_t index = cmatched->GetHitIndex(ii);
 583                 AliESDTOFHit *hitOld = (AliESDTOFHit *) esdTOFHitArr->At(index);
 584                 Int_t index_new = fHitsESD->GetEntriesFast();
 585                 AliESDTOFHit *hitNew = new( (*fHitsESD)[index_new] ) AliESDTOFHit(*hitOld);
 586                 hitNew->SetESDTOFClusterIndex(currentpos);
 587                 clnew->SetHitIndex(ii,index_new);
 588               }
 589
 590               fWrittenInPos[clind[i]] = currentpos;
 591               cmatched = clnew; // update the new one added to the ESDEvent
 592             }
 593
 594             if(cmatched->GetNMatchableTracks() < AliESDTOFCluster::kMaxMatches){
 595               cmatched->Update(t->GetID(),dist3d[0],dist3d[1],dist3d[2],fTrackPos[3][istep],timesCurrent);//x,y,z -> tracking RF
 596               t->AddTOFcluster(currentpos);
 597               t->SetStatus(AliESDtrack::kTOFout);
 598             }
 599           }
 600           AliDebug(2,Form(" dist3dLoc[0] = %f, dist3dLoc[1] = %f, dist3dLoc[2] = %f ",dist3d[0],dist3d[1],dist3d[2]));
 601
 602           nfound++;
 603
 604           // ***** NEW *****
 605         }//end if accept
 606
 607       } //end for on the clusters
 608     } //end for on the steps
 609
 610
 611     if (nfound == 0 ) {
 612       AliDebug(1,Form(" No matchable track points for the track number %d",iseed));
 613       fnunmatch++;
 614       continue;
 615     }
 616
 617     AliDebug(1,Form(" Number of track points for the track number %d: %d",iseed,nfound));
 618
 619     Int_t nMatchedClusters = t->GetNTOFclusters();
 620
 621     if (nMatchedClusters==0) {
 622       AliDebug(1,Form("Reconstructed track %d doesn't match any TOF cluster", iseed));
 623       fnunmatch++;
 624       continue;
 625     }
 626
 627     AliDebug(1,Form(" %d - matched (%d)",iseed,nMatchedClusters));
 628
 629     fnmatch++;
 630
 631     /*
 632     AliTOFcluster cTOF = AliTOFcluster(volIdClus,
 633     (Float_t)posClus[0],(Float_t)posClus[1],(Float_t)posClus[2],
 634     (Float_t)covClus[0],(Float_t)covClus[1],(Float_t)covClus[2],
 635     (Float_t)covClus[3],(Float_t)covClus[4],(Float_t)covClus[5],
 636     tofLabels,volIndices,parClu,kTRUE,index[i]);
 637
 638     // Fill the track residual histograms.
 639     FillResiduals(trackTOFin,c,kFALSE);
 640     */
 641   } // loop on fSeeds
 642
 643 }
 644 //_________________________________________________________________________
 645 Int_t AliTOFtrackerV2::LoadClusters(TTree *cTree) {
 646   //--------------------------------------------------------------------
 647   //This function loads the TOF clusters
 648   //--------------------------------------------------------------------
 649
 650   TBranch *branch=cTree->GetBranch("TOF");
 651   if (!branch) {
 652     AliError("can't get the branch with the TOF clusters !");
 653     return 1;
 654   }
 655
 656   static TClonesArray dummy("AliTOFcluster",10000);
 657   dummy.Clear();
 658   TClonesArray *clusters=&dummy;
 659   branch->SetAddress(&clusters);
 660
 661   cTree->GetEvent(0);
 662   Int_t ncl =clusters->GetEntriesFast();
 663   AliInfo(Form("Number of clusters: %d",ncl));
 664
 665   fN = ncl; // set cluster counter
 666
 667   for(Int_t i=0;i < ncl;i++) // reset position of clusters in ESD
 668     fWrittenInPos[i] = -1;
 669
 670   if(ncl==0){
 671     return 0;
 672   }
 673
 674   for (Int_t i=0; i<ncl; i++) {
 675     AliTOFcluster *c=(AliTOFcluster*)clusters->UncheckedAt(i);
 676
 677     /*
 678     Int_t ind[5];
 679     ind[0]=c->GetDetInd(0);
 680     ind[1]=c->GetDetInd(1);
 681     ind[2]=c->GetDetInd(2);
 682     ind[3]=c->GetDetInd(3);
 683     ind[4]=c->GetDetInd(4);
 684     Int_t calindex = AliTOFGeometry::GetIndex(ind);
 685     Int_t tofLabels[3]={c->GetLabel(0),c->GetLabel(1),c->GetLabel(2)};
 686     */
 687
 688     fClusters[i] = c;
 689   }
 690
 691   return 0;
 692 }
 693 //_________________________________________________________________________
 694 void AliTOFtrackerV2::UnloadClusters() {
 695   //--------------------------------------------------------------------
 696   //This function unloads TOF clusters
 697   //--------------------------------------------------------------------
 698
 699   // don't delete TOF clusters here because they should be written
 700 }
 701
 702 //_________________________________________________________________________
 703 Int_t AliTOFtrackerV2::FindClusterIndex(Double_t z) const {
 704   //--------------------------------------------------------------------
 705   // This function returns the index of the nearest cluster
 706   //--------------------------------------------------------------------
 707   TClonesArray* TOFClArr = fClustersESD;; // use temporary array
 708   Int_t n = TOFClArr->GetEntriesFast();
 709
 710   if (n==0) return 0;
 711   if (z <= ((AliESDTOFCluster *) TOFClArr->At(0))->GetZ()) return 0;
 712   if (z > ((AliESDTOFCluster *) TOFClArr->At(n-1))->GetZ()) return n;
 713   Int_t b=0, e=n-1, m=(b+e)/2;
 714   for (; b<e; m=(b+e)/2) {
 715     if (z > ((AliESDTOFCluster *) TOFClArr->At(m))->GetZ()) b=m+1;
 716
 717     else e=m;
 718   }
 719   return m;
 720 }
 721 //_________________________________________________________________________
 722 Bool_t AliTOFtrackerV2::GetTrackPoint(Int_t index, AliTrackPoint& p) const
 723 {
 724   // Get track space point with index i
 725   // Coordinates are in the global system
 726   TClonesArray* esdTOFClArr = fEvent->GetESDTOFClusters();
 727   AliESDTOFCluster *cl = (AliESDTOFCluster *) esdTOFClArr->At(index);
 728
 729   Float_t xyz[3];
 730   xyz[0] = cl->GetR()*TMath::Cos(cl->GetPhi());
 731   xyz[1] = cl->GetR()*TMath::Sin(cl->GetPhi());
 732   xyz[2] = cl->GetZ();
 733   Float_t phiangle = (Int_t(cl->GetPhi()*TMath::RadToDeg()/20.)+0.5)*20.*TMath::DegToRad();
 734   Float_t sinphi = TMath::Sin(phiangle), cosphi = TMath::Cos(phiangle);
 735   Int_t tofChannel=cl->GetTOFchannel();
 736   Int_t ind[5]; AliTOFGeometry::GetVolumeIndices(tofChannel,ind);
 737   Float_t tiltangle = AliTOFGeometry::GetAngles(ind[1],ind[2])*TMath::DegToRad();
 738   Float_t sinth = TMath::Sin(tiltangle), costh = TMath::Cos(tiltangle);
 739   Float_t sigmay2 = AliTOFGeometry::XPad()*AliTOFGeometry::XPad()/12.;
 740   Float_t sigmaz2 = AliTOFGeometry::ZPad()*AliTOFGeometry::ZPad()/12.;
 741   Float_t cov[6];
 742   cov[0] = sinphi*sinphi*sigmay2 + cosphi*cosphi*sinth*sinth*sigmaz2;
 743   cov[1] = -sinphi*cosphi*sigmay2 + sinphi*cosphi*sinth*sinth*sigmaz2;
 744   cov[2] = -cosphi*sinth*costh*sigmaz2;
 745   cov[3] = cosphi*cosphi*sigmay2 + sinphi*sinphi*sinth*sinth*sigmaz2;
 746   cov[4] = -sinphi*sinth*costh*sigmaz2;
 747   cov[5] = costh*costh*sigmaz2;
 748   p.SetXYZ(xyz[0],xyz[1],xyz[2],cov);
 749
 750   // Detector numbering scheme
 751   Int_t nSector = AliTOFGeometry::NSectors();
 752   Int_t nPlate  = AliTOFGeometry::NPlates();
 753   Int_t nStripA = AliTOFGeometry::NStripA();
 754   Int_t nStripB = AliTOFGeometry::NStripB();
 755   Int_t nStripC = AliTOFGeometry::NStripC();
 756
 757   Int_t isector = ind[0];//cl->GetDetInd(0);
 758   if (isector >= nSector)
 759     AliError(Form("Wrong sector number in TOF (%d) !",isector));
 760   Int_t iplate = ind[1];//cl->GetDetInd(1);
 761   if (iplate >= nPlate)
 762     AliError(Form("Wrong plate number in TOF (%d) !",iplate));
 763   Int_t istrip = ind[2];//cl->GetDetInd(2);
 764
 765   Int_t stripOffset = 0;
 766   switch (iplate) {
 767   case 0:
 768     stripOffset = 0;
 769     break;
 770   case 1:
 771     stripOffset = nStripC;
 772     break;
 773   case 2:
 774     stripOffset = nStripC+nStripB;
 775     break;
 776   case 3:
 777     stripOffset = nStripC+nStripB+nStripA;
 778     break;
 779   case 4:
 780     stripOffset = nStripC+nStripB+nStripA+nStripB;
 781     break;
 782   default:
 783     AliError(Form("Wrong plate number in TOF (%d) !",iplate));
 784     break;
 785   };
 786
 787   Int_t idet = (2*(nStripC+nStripB)+nStripA)*isector +
 788                stripOffset +
 789                istrip;
 790   UShort_t volid = AliGeomManager::LayerToVolUID(AliGeomManager::kTOF,idet);
 791   p.SetVolumeID((UShort_t)volid);
 792   return kTRUE;
 793 }
 794 //_________________________________________________________________________
 795
 796 Float_t AliTOFtrackerV2::CorrectTimeWalk( Float_t dist, Float_t tof) const {
 797
 798   //dummy, for the moment
 799   Float_t tofcorr=0.;
 800   if(dist<AliTOFGeometry::ZPad()*0.5){
 801     tofcorr=tof;
 802     //place here the actual correction
 803   }else{
 804     tofcorr=tof;
 805   }
 806   return tofcorr;
 807 }
 808 //_________________________________________________________________________
 809 void AliTOFtrackerV2::Clusterize(){
 810   Int_t detId[5];
 811
 812   // Load 1 hit in 1 cluster (ESD)
 813   TClonesArray* TOFClArr = fClustersESD;// use a temporary copy //fEvent->GetESDTOFClusters();
 814   TClonesArray* esdTOFHitArr = fEvent->GetESDTOFHits();
 815
 816   if(TOFClArr->GetEntriesFast()) TOFClArr->Clear();
 817   if(esdTOFHitArr->GetEntriesFast()) esdTOFHitArr->Clear();
 818
 819   AliESDTOFCluster *c1 = NULL;
 820   AliESDTOFCluster *c2 = NULL;
 821
 822   for(Int_t i=0; i < fN;i++){
 823     AliTOFcluster *c = fClusters[i];
 824     Int_t ind[5];
 825     ind[0]=c->GetDetInd(0);
 826     ind[1]=c->GetDetInd(1);
 827     ind[2]=c->GetDetInd(2);
 828     ind[3]=c->GetDetInd(3);
 829     ind[4]=c->GetDetInd(4);
 830     Int_t calindex = AliTOFGeometry::GetIndex(ind);
 831     Int_t tofLabels[3]={c->GetLabel(0),c->GetLabel(1),c->GetLabel(2)};
 832
 833     new ( (*esdTOFHitArr)[i] ) AliESDTOFHit( AliTOFGeometry::TdcBinWidth()*c->GetTDC(),
 834                               AliTOFGeometry::TdcBinWidth()*c->GetTDCRAW(),
 835                               AliTOFGeometry::ToTBinWidth()*c->GetToT()*1E-3,
 836                               calindex,tofLabels,c->GetL0L1Latency(),
 837                               c->GetDeltaBC(),i,c->GetZ(),c->GetR(),c->GetPhi() );
 838
 839     c1 =  new( (*TOFClArr)[i] ) AliESDTOFCluster(i);
 840     c1->SetEvent(fEvent);
 841     c1->SetStatus( c->GetStatus() );
 842     c1->SetESDID(i);
 843     //
 844     // register hits in the cluster
 845     c1->AddESDTOFHitIndex(i);
 846
 847   }
 848   // start to merge clusters
 849   Int_t chan1,chan2,chan3;
 850   Int_t strip1,strip2;
 851   Int_t iphi,iphi2,iphi3;
 852   Int_t ieta,ieta2,ieta3;
 853
 854   for(Int_t i=0; i <  TOFClArr->GetEntriesFast()-1;i++){
 855     c1=(AliESDTOFCluster *) TOFClArr->At(i);
 856     if(!c1->GetStatus()) continue;
 857
 858     chan1 = c1->GetTOFchannel(0);
 859     AliTOFGeometry::GetVolumeIndices(chan1, detId); // Get volume index from channel index
 860
 861     ieta = detId[2]/*strip*/*2 + detId[3]/*pad Z*/;
 862     if(detId[1]/*module*/ == 0) ieta += 0;
 863     else if(detId[1] == 1) ieta += 38;
 864     else if(detId[1] == 2) ieta += 76;
 865     else if(detId[1] == 3) ieta += 106;
 866     else if(detId[1] == 4) ieta += 144;
 867     iphi = detId[0]/*phi sector*/*48 + detId[4]/*pad x*/;
 868
 869     chan2=chan1;
 870     if(c1->GetNTOFhits() > 1){
 871       chan2 = c1->GetTOFchannel(1);
 872       AliTOFGeometry::GetVolumeIndices(chan2, detId); // Get volume index from channel index
 873
 874       ieta2 = detId[2]/*strip*/*2 + detId[3]/*pad Z*/;
 875       if(detId[1]/*module*/ == 0) ieta2 += 0;
 876       else if(detId[1] == 1) ieta2 += 38;
 877       else if(detId[1] == 2) ieta2 += 76;
 878       else if(detId[1] == 3) ieta2 += 106;
 879       else if(detId[1] == 4) ieta2 += 144;
 880       iphi2 = detId[0]/*phi sector*/*48 + detId[4]/*pad x*/;
 881     }
 882     else{
 883       iphi2=iphi;
 884       ieta2=ieta;
 885     }
 886
 887     // 1 and 2 belong now to the first cluster, 3 to the second one
 888
 889     strip1 = chan1/96;
 890     for(Int_t j=i+1; j < TOFClArr->GetEntriesFast();j++){
 891       c2=(AliESDTOFCluster *) TOFClArr->At(j);
 892       if(!c2->GetStatus()) continue;
 893
 894       chan3 = c2->GetTOFchannel();
 895
 896       // check if the two TOF hits are in the same strip
 897       strip2 = chan3/96;
 898       if(strip1 != strip2) continue;
 899
 900       AliTOFGeometry::GetVolumeIndices(chan3, detId); // Get volume index from channel index
 901       ieta3 = detId[2]/*strip*/*2 + detId[3]/*pad Z*/;
 902       if(detId[1]/*module*/ == 0) ieta3 += 0;
 903       else if(detId[1] == 1) ieta3 += 38;
 904       else if(detId[1] == 2) ieta3 += 76;
 905       else if(detId[1] == 3) ieta3 += 106;
 906       else if(detId[1] == 4) ieta3 += 144;
 907       iphi3 = detId[0]/*phi sector*/*48 + detId[4]/*pad x*/;
 908
 909
 910       if(ieta3-ieta > 2) j = fN; // because cluster are order along Z so you can skip all the rest, go to the next one ("i+1")
 911
 912       // check if the fired pad are close in space
 913       if((TMath::Abs(iphi-iphi3)>1 && TMath::Abs(iphi2-iphi3)>1) || (TMath::Abs(ieta-ieta3)>1 && TMath::Abs(ieta2-ieta3)>1))
 914         continue; // double checks
 915
 916       // check if the TOF time are close enough to be merged
 917       if(TMath::Abs(c1->GetTime() - c2->GetTime()) > 500/*in ps*/) continue;
 918
 919       // merge them
 920       MergeClusters(i,j);
 921
 922       // new hit is added as a second hit for the first cluster
 923       iphi2 = iphi3;
 924       ieta2 = ieta3;
 925     }
 926   }
 927 }
 928
 929 void AliTOFtrackerV2::MergeClusters(Int_t i,Int_t j){
 930   TClonesArray* TOFClArr = fClustersESD;// use a temporary copy //fEvent->GetESDTOFClusters();
 931
 932   if(i == j){
 933     AliInfo("No TOF cluster mergine possible (cannot merge a cluster with itself)");
 934     return;
 935   }
 936
 937   if(i > j){ // check right order
 938     Int_t k=i;
 939     i=j;
 940     j=k;
 941   }
 942
 943   Int_t last = TOFClArr->GetEntriesFast()-1;
 944
 945   if(j > last){
 946     AliInfo("No TOF cluster mergine possible (cluster not available)");
 947     return;
 948   }
 949
 950   AliESDTOFCluster *c1 = (AliESDTOFCluster *) TOFClArr->At(i);
 951   AliESDTOFCluster *c2 = (AliESDTOFCluster *) TOFClArr->At(j);
 952
 953   if(c2->GetNMatchableTracks()){
 954     AliInfo("No TOF cluster mergine possible (cluster already matched)");
 955     return; // cannot merge a cluster already matched
 956   }
 957
 958   Int_t nhit1 = c1->GetNTOFhits();
 959   Int_t nhit2 = c2->GetNTOFhits();
 960
 961   if(nhit1+nhit2 >= AliESDTOFCluster::kMaxHits)
 962     {
 963       AliInfo("No TOF cluster mergine possible (too many hits)");
 964       return;
 965     }
 966
 967   for(Int_t k=0;k < nhit2;k++){// add hits in c2 to c1
 968     c1->AddESDTOFHitIndex(c2->GetHitIndex(k));
 969
 970     // ID re-setting for hits not needed (done when matching is found)
 971   }
 972
 973   // remove c2 from array
 974   if(j == last) delete TOFClArr->RemoveAt(j);
 975   else{
 976     for(Int_t ii=j;ii < last;ii++){
 977       AliESDTOFCluster *old= (AliESDTOFCluster *) TOFClArr->At(ii);
 978       if (!old) {AliFatal(Form("NULL pointer for TOF cluster %d",ii));}
 979       AliESDTOFCluster *replace= (AliESDTOFCluster *) TOFClArr->At(ii+1);
 980       if (!replace) {AliFatal(Form("NULL pointer for TOF cluster %d",ii+1));}
 981       *old = *replace;
 982       old->SetESDID(j);
 983     }
 984     delete TOFClArr->RemoveAt(last);
 985   }
 986
 987 }