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