TOF/AliTOFtracker.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 // AliTOFtracker 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 <TSeqCollection.h>
  33 #include <TClonesArray.h>
  34 #include <TObjArray.h>
  35 #include <TGeoManager.h>
  36 #include <TTree.h>
  37 #include <TFile.h>
  38 #include <TH2F.h>
  39
  40 #include "AliGeomManager.h"
  41 #include "AliESDtrack.h"
  42 #include "AliESDEvent.h"
  43 #include "AliESDpid.h"
  44 #include "AliLog.h"
  45 #include "AliTrackPointArray.h"
  46 #include "AliCDBManager.h"
  47
  48 //#include "AliTOFpidESD.h"
  49 #include "AliTOFRecoParam.h"
  50 #include "AliTOFReconstructor.h"
  51 #include "AliTOFcluster.h"
  52 #include "AliTOFGeometry.h"
  53 #include "AliTOFtracker.h"
  54 #include "AliTOFtrack.h"
  55
  56 extern TGeoManager *gGeoManager;
  57 //extern TROOT *gROOT;
  58
  59
  60 ClassImp(AliTOFtracker)
  61
  62 //_____________________________________________________________________________
  63 AliTOFtracker::AliTOFtracker():
  64   fkRecoParam(0x0),
  65   fGeom(0x0),
  66   fN(0),
  67   fNseeds(0),
  68   fNseedsTOF(0),
  69   fngoodmatch(0),
  70   fnbadmatch(0),
  71   fnunmatch(0),
  72   fnmatch(0),
  73   fESDEv(0),
  74   fTracks(new TClonesArray("AliTOFtrack")),
  75   fSeeds(new TObjArray(100)),
  76   fTOFtrackPoints(new TObjArray(10)),
  77   fHDigClusMap(0x0),
  78   fHDigNClus(0x0),
  79   fHDigClusTime(0x0),
  80   fHDigClusToT(0x0),
  81   fHRecNClus(0x0),
  82   fHRecDist(0x0),
  83   fHRecSigYVsP(0x0),
  84   fHRecSigZVsP(0x0),
  85   fHRecSigYVsPWin(0x0),
  86   fHRecSigZVsPWin(0x0),
  87   fCalTree(0x0),
  88   fIch(-1),
  89   fToT(-1.),
  90   fTime(-1.),
  91   fExpTimePi(-1.),
  92   fExpTimeKa(-1.),
  93   fExpTimePr(-1.),
  94   fNTOFmatched(0)
  95 {
  96   //AliTOFtracker main Ctor
  97
  98   for (Int_t ii=0; ii<kMaxCluster; ii++) fClusters[ii]=0x0;
  99
 100   // Getting the geometry
 101   fGeom = new AliTOFGeometry();
 102   /* RS?
 103   for(Int_t i=0; i< 20000;i++){
 104     fClusterESD[i] = NULL;
 105     fHit[i] = NULL;
 106   }
 107   */
 108   InitCheckHists();
 109
 110 }
 111 //_____________________________________________________________________________
 112 AliTOFtracker::~AliTOFtracker() {
 113   //
 114   // Dtor
 115   //
 116
 117   SaveCheckHists();
 118
 119   if(!(AliCDBManager::Instance()->GetCacheFlag())){
 120     delete fkRecoParam;
 121   }
 122   delete fGeom;
 123   delete fHDigClusMap;
 124   delete fHDigNClus;
 125   delete fHDigClusTime;
 126   delete fHDigClusToT;
 127   delete fHRecNClus;
 128   delete fHRecDist;
 129   delete fHRecSigYVsP;
 130   delete fHRecSigZVsP;
 131   delete fHRecSigYVsPWin;
 132   delete fHRecSigZVsPWin;
 133   delete fCalTree;
 134   if (fTracks){
 135     fTracks->Delete();
 136     delete fTracks;
 137     fTracks=0x0;
 138   }
 139   if (fSeeds){
 140     fSeeds->Delete();
 141     delete fSeeds;
 142     fSeeds=0x0;
 143   }
 144   if (fTOFtrackPoints){
 145     fTOFtrackPoints->Delete();
 146     delete fTOFtrackPoints;
 147     fTOFtrackPoints=0x0;
 148   }
 149
 150   for (Int_t ii=0; ii<kMaxCluster; ii++)
 151     if (fClusters[ii]) fClusters[ii]->Delete();
 152
 153   /* RS?
 154      for(Int_t i=0; i< 20000;i++){
 155      if(fClusterESD[i]){
 156      delete fClusterESD[i];
 157      fClusterESD[i] = NULL;
 158      }
 159      if(fHit[i]){
 160      delete fHit[i];
 161      fHit[i] = NULL;
 162      }
 163      }
 164   */
 165
 166 }
 167 //_____________________________________________________________________________
 168 void AliTOFtracker::GetPidSettings(AliESDpid *esdPID) {
 169   //
 170   // Sets TOF resolution from RecoParams
 171   //
 172   if (fkRecoParam)
 173     esdPID->GetTOFResponse().SetTimeResolution(fkRecoParam->GetTimeResolution());
 174   else
 175     AliWarning("fkRecoParam not yet set; cannot set PID settings");
 176 }
 177 //_____________________________________________________________________________
 178 Int_t AliTOFtracker::PropagateBack(AliESDEvent * const event) {
 179   //
 180   // Gets seeds from ESD event and Match with TOF Clusters
 181   //
 182   fESDEv = event;
 183   //
 184   if (fN==0) {
 185     AliInfo("No TOF recPoints to be matched with reconstructed tracks");
 186     return 0;
 187   }
 188
 189   // initialize RecoParam for current event
 190   AliDebug(1,"Initializing params for TOF");
 191
 192   fkRecoParam = AliTOFReconstructor::GetRecoParam();  // instantiate reco param from STEER...
 193
 194   if (fkRecoParam == 0x0) {
 195     AliFatal("No Reco Param found for TOF!!!");
 196   }
 197   //fkRecoParam->Dump();
 198   //if(fkRecoParam->GetApplyPbPbCuts())fkRecoParam=fkRecoParam->GetPbPbparam();
 199   //fkRecoParam->PrintParameters();
 200
 201   /* RS?
 202   // create clusters from hit
 203   for(Int_t i=0;i < fNTOFmatched;i++){
 204     fClusterESD[i] = new AliESDTOFCluster(fHit[i],event);
 205     fClusterESD[i]->SetStatus(fClusters[i]->GetStatus());
 206   }
 207   */
 208   //Initialise some counters
 209
 210   fNseeds=0;
 211   fNseedsTOF=0;
 212   fngoodmatch=0;
 213   fnbadmatch=0;
 214   fnunmatch=0;
 215   fnmatch=0;
 216
 217   Int_t ntrk=event->GetNumberOfTracks();
 218   fNseeds = ntrk;
 219
 220
 221   //Load ESD tracks into a local Array of ESD Seeds
 222   for (Int_t i=0; i<fNseeds; i++){
 223     fSeeds->AddLast(event->GetTrack(i));
 224     //    event->GetTrack(i)->SetESDEvent(event); // RS: Why this is needed? The event is already set
 225   }
 226   //Prepare ESD tracks candidates for TOF Matching
 227   CollectESD();
 228
 229   if (fNseeds==0 || fNseedsTOF==0) {
 230     AliInfo("No seeds to try TOF match");
 231     fSeeds->Clear();
 232     fTracks->Clear();
 233     return 0 ;
 234   }
 235
 236   //First Step with Strict Matching Criterion
 237   MatchTracks(0);
 238
 239   //Second Step with Looser Matching Criterion
 240   MatchTracks(1);
 241
 242   //Third Step without kTOFout flag (just to update clusters)
 243   MatchTracks(2);
 244
 245   //RS?  event->SetTOFcluster(fNTOFmatched,fClusterESD);
 246
 247   if (fN==0) {
 248     AliInfo("No TOF recPoints to be matched with reconstructed tracks");
 249     fSeeds->Clear();
 250     fTracks->Clear();
 251     return 0;
 252   }
 253
 254   AliInfo(Form("Number of matched tracks = %d (good = %d, bad = %d)",fnmatch,fngoodmatch,fnbadmatch));
 255
 256   //Update the matched ESD tracks
 257
 258   for (Int_t i=0; i<ntrk; i++) {
 259     // RS: This is a bogus code since t and seed are the same object
 260     //    AliESDtrack *t=event->GetTrack(i);
 261     AliESDtrack *seed =(AliESDtrack*)fSeeds->At(i);
 262     if ( (seed->GetStatus()&AliESDtrack::kTOFin)!=0 ) {
 263       //t->SetStatus(AliESDtrack::kTOFin);
 264       //if(seed->GetTOFsignal()>0){
 265       if ( (seed->GetStatus()&AliESDtrack::kTOFout)!=0 ) {
 266         //t->SetStatus(AliESDtrack::kTOFout);
 267         //t->SetTOFsignal(seed->GetTOFsignal());
 268         //t->SetTOFcluster(seed->GetTOFcluster());
 269         //t->SetTOFsignalToT(seed->GetTOFsignalToT());
 270         //t->SetTOFsignalRaw(seed->GetTOFsignalRaw());
 271         //t->SetTOFsignalDz(seed->GetTOFsignalDz());
 272         //t->SetTOFsignalDx(seed->GetTOFsignalDx());
 273         //t->SetTOFDeltaBC(seed->GetTOFDeltaBC());
 274         //t->SetTOFL0L1(seed->GetTOFL0L1());
 275         //t->SetTOFCalChannel(seed->GetTOFCalChannel());
 276         Int_t tlab[3]; seed->GetTOFLabel(tlab);
 277         //t->SetTOFLabel(tlab);
 278
 279         Double_t alphaA = (Double_t)seed->GetAlpha();
 280         Double_t xA = (Double_t)seed->GetX();
 281         Double_t yA = (Double_t)seed->GetY();
 282         Double_t zA = (Double_t)seed->GetZ();
 283         Double_t p1A = (Double_t)seed->GetSnp();
 284         Double_t p2A = (Double_t)seed->GetTgl();
 285         Double_t p3A = (Double_t)seed->GetSigned1Pt();
 286         const Double_t *covA = (Double_t*)seed->GetCovariance();
 287
 288         // Make attention, please:
 289         //      AliESDtrack::fTOFInfo array does not be stored in the AliESDs.root file
 290         //      it is there only for a check during the reconstruction step.
 291         Float_t info[10]; seed->GetTOFInfo(info);
 292         //t->SetTOFInfo(info);
 293         AliDebug(3,Form(" distance=%f; residual in the pad reference frame: dX=%f, dZ=%f", info[0],info[1],info[2]));
 294
 295         // Check done:
 296         //       by calling the AliESDtrack::UpdateTrackParams,
 297         //       the current track parameters are changed
 298         //       and it could cause refit problems.
 299         //       We need to update only the following track parameters:
 300         //            the track length and expected times.
 301         //       Removed AliESDtrack::UpdateTrackParams call
 302         //       Called AliESDtrack::SetIntegratedTimes(...) and
 303         //       AliESDtrack::SetIntegratedLength() routines.
 304         /*
 305           AliTOFtrack *track = new AliTOFtrack(*seed);
 306           t->UpdateTrackParams(track,AliESDtrack::kTOFout); // to be checked - AdC
 307           delete track;
 308           Double_t time[AliPID::kSPECIESC]; t->GetIntegratedTimes(time);
 309         */
 310
 311         Double_t time[AliPID::kSPECIESC]; seed->GetIntegratedTimes(time,AliPID::kSPECIESC);
 312         //t->SetIntegratedTimes(time);
 313
 314         //Double_t length =  seed->GetIntegratedLength();
 315         //t->SetIntegratedLength(length);
 316
 317         Double_t alphaB = (Double_t)seed->GetAlpha();
 318         Double_t xB = (Double_t)seed->GetX();
 319         Double_t yB = (Double_t)seed->GetY();
 320         Double_t zB = (Double_t)seed->GetZ();
 321         Double_t p1B = (Double_t)seed->GetSnp();
 322         Double_t p2B = (Double_t)seed->GetTgl();
 323         Double_t p3B = (Double_t)seed->GetSigned1Pt();
 324         const Double_t *covB = (Double_t*)seed->GetCovariance();
 325         AliDebug(3,"Track params -now(before)-:");
 326         AliDebug(3,Form("    X: %f(%f), Y: %f(%f), Z: %f(%f) --- alpha: %f(%f)",
 327                         xB,xA,
 328                         yB,yA,
 329                         zB,zA,
 330                         alphaB,alphaA));
 331         AliDebug(3,Form("    p1: %f(%f), p2: %f(%f), p3: %f(%f)",
 332                         p1B,p1A,
 333                         p2B,p2A,
 334                         p3B,p3A));
 335         AliDebug(3,Form("    cov1: %f(%f), cov2: %f(%f), cov3: %f(%f)"
 336                         " cov4: %f(%f), cov5: %f(%f), cov6: %f(%f)"
 337                         " cov7: %f(%f), cov8: %f(%f), cov9: %f(%f)"
 338                         " cov10: %f(%f), cov11: %f(%f), cov12: %f(%f)"
 339                         " cov13: %f(%f), cov14: %f(%f), cov15: %f(%f)",
 340                         covB[0],covA[0],
 341                         covB[1],covA[1],
 342                         covB[2],covA[2],
 343                         covB[3],covA[3],
 344                         covB[4],covA[4],
 345                         covB[5],covA[5],
 346                         covB[6],covA[6],
 347                         covB[7],covA[7],
 348                         covB[8],covA[8],
 349                         covB[9],covA[9],
 350                         covB[10],covA[10],
 351                         covB[11],covA[11],
 352                         covB[12],covA[12],
 353                         covB[13],covA[13],
 354                         covB[14],covA[14]
 355                         ));
 356         AliDebug(2,Form(" TOF params: %6d  %f %f %f %f %f %6d %3d  %f",
 357                         i,
 358                         seed->GetTOFsignalRaw(),
 359                         seed->GetTOFsignal(),
 360                         seed->GetTOFsignalToT(),
 361                         seed->GetTOFsignalDz(),
 362                         seed->GetTOFsignalDx(),
 363                         seed->GetTOFCalChannel(),
 364                         seed->GetTOFcluster(),
 365                         seed->GetIntegratedLength()));
 366         AliDebug(2,Form(" %f %f %f %f %f %f %f %f %f",
 367                         time[0], time[1], time[2], time[3], time[4], time[5], time[6], time[7], time[8]));
 368       }
 369     }
 370   }
 371   /* RS?
 372   if(fNTOFmatched){
 373     Int_t *matchmap = new Int_t[fNTOFmatched];
 374     event->SetTOFcluster(fNTOFmatched,fClusterESD,matchmap);
 375     for (Int_t i=0; i<ntrk; i++) { // remapping after TOF matching selection
 376       AliESDtrack *t=event->GetTrack(i);
 377       t->ReMapTOFcluster(fNTOFmatched,matchmap);
 378     }
 379
 380     delete[] matchmap;
 381   }
 382   */
 383
 384   //Make TOF PID
 385   // Now done in AliESDpid
 386   // fPid->MakePID(event,timeZero);
 387
 388   fSeeds->Clear();
 389   //fTracks->Delete();
 390   fTracks->Clear();
 391   return 0;
 392
 393 }
 394 //_________________________________________________________________________
 395 void AliTOFtracker::CollectESD() {
 396    //prepare the set of ESD tracks to be matched to clusters in TOF
 397
 398   Int_t seedsTOF1=0;
 399   Int_t seedsTOF3=0;
 400   Int_t seedsTOF2=0;
 401
 402   TClonesArray &aTOFTrack = *fTracks;
 403   for (Int_t i=0; i<fNseeds; i++) {
 404
 405     AliESDtrack *t =(AliESDtrack*)fSeeds->At(i);
 406     if ((t->GetStatus()&AliESDtrack::kTPCout)==0)continue;
 407
 408     AliTOFtrack *track = new AliTOFtrack(*t); // New
 409     Float_t x = (Float_t)track->GetX(); //New
 410
 411     // TRD 'good' tracks
 412     if ( ( (t->GetStatus()&AliESDtrack::kTRDout)!=0 ) ) {
 413
 414       AliDebug(1,Form(" Before propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track->GetIntegratedLength()));
 415
 416       // TRD 'good' tracks, already propagated at 371 cm
 417       if( x >= AliTOFGeometry::Rmin() ) {
 418
 419         if  ( track->PropagateToInnerTOF() ) {
 420
 421           AliDebug(1,Form(" TRD propagated track till rho = %fcm."
 422                           " And then the track has been propagated till rho = %fcm.",
 423                           x, (Float_t)track->GetX()));
 424
 425           track->SetSeedIndex(i);
 426           t->UpdateTrackParams(track,AliESDtrack::kTOFin);
 427           new(aTOFTrack[fNseedsTOF]) AliTOFtrack(*track);
 428           fNseedsTOF++;
 429           seedsTOF1++;
 430
 431           AliDebug(1,Form(" After propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track->GetIntegratedLength()));
 432         }
 433         delete track;
 434
 435       }
 436       else { // TRD 'good' tracks, propagated rho<371cm
 437
 438         if  ( track->PropagateToInnerTOF() ) {
 439
 440           AliDebug(1,Form(" TRD propagated track till rho = %fcm."
 441                           " And then the track has been propagated till rho = %fcm.",
 442                           x, (Float_t)track->GetX()));
 443
 444           track->SetSeedIndex(i);
 445           t->UpdateTrackParams(track,AliESDtrack::kTOFin);
 446           new(aTOFTrack[fNseedsTOF]) AliTOFtrack(*track);
 447           fNseedsTOF++;
 448           seedsTOF3++;
 449
 450           AliDebug(1,Form(" After propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track->GetIntegratedLength()));
 451         }
 452         delete track;
 453
 454       }
 455       //delete track;
 456     }
 457
 458     else { // Propagate the rest of TPCbp
 459
 460       AliDebug(1,Form(" Before propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track->GetIntegratedLength()));
 461
 462       if ( track->PropagateToInnerTOF() ) {
 463
 464         AliDebug(1,Form(" TPC propagated track till rho = %fcm."
 465                         " And then the track has been propagated till rho = %fcm.",
 466                         x, (Float_t)track->GetX()));
 467
 468         track->SetSeedIndex(i);
 469         t->UpdateTrackParams(track,AliESDtrack::kTOFin);
 470         new(aTOFTrack[fNseedsTOF]) AliTOFtrack(*track);
 471         fNseedsTOF++;
 472         seedsTOF2++;
 473
 474         AliDebug(1,Form(" After propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track->GetIntegratedLength()));
 475       }
 476       delete track;
 477     }
 478   }
 479
 480   AliInfo(Form("Number of TOF seeds = %d (kTRDout371 = %d, kTRDoutLess371 = %d, !kTRDout = %d)",fNseedsTOF,seedsTOF1,seedsTOF3,seedsTOF2));
 481
 482   // Sort according uncertainties on track position
 483   fTracks->Sort();
 484
 485 }
 486
 487 //_________________________________________________________________________
 488 void AliTOFtracker::MatchTracks( Int_t mLastStep){
 489
 490   // Parameters used/regulating the reconstruction
 491
 492   //static Float_t corrLen=0.;//0.75;
 493   static Float_t detDepth=18.;
 494   static Float_t padDepth=0.5;
 495
 496   const Float_t kSpeedOfLight= 2.99792458e-2; // speed of light [cm/ps]
 497   const Float_t kTimeOffset = 0.; // time offset for tracking algorithm [ps]
 498
 499   Float_t dY=AliTOFGeometry::XPad();
 500   Float_t dZ=AliTOFGeometry::ZPad();
 501
 502   Float_t sensRadius = fkRecoParam->GetSensRadius();
 503   Float_t stepSize   = fkRecoParam->GetStepSize();
 504   Float_t scaleFact  = fkRecoParam->GetWindowScaleFact();
 505   Float_t dyMax=fkRecoParam->GetWindowSizeMaxY();
 506   Float_t dzMax=fkRecoParam->GetWindowSizeMaxZ();
 507   Float_t dCut=fkRecoParam->GetDistanceCut();
 508   if (dCut==3. && fNseedsTOF<=10) {
 509     dCut=10.;
 510     AliInfo(Form("Matching window=%f, since low multiplicity event (fNseedsTOF=%d)",
 511                  dCut, fNseedsTOF));
 512   }
 513   if(mLastStep == 2)
 514     dCut=10.;
 515
 516
 517   Double_t maxChi2=fkRecoParam->GetMaxChi2TRD();
 518   Bool_t timeWalkCorr    = fkRecoParam->GetTimeWalkCorr();
 519   if(!mLastStep){
 520     AliDebug(1,"++++++++++++++TOF Reconstruction Parameters:++++++++++++");
 521     AliDebug(1,Form("TOF sens radius: %f",sensRadius));
 522     AliDebug(1,Form("TOF step size: %f",stepSize));
 523     AliDebug(1,Form("TOF Window scale factor: %f",scaleFact));
 524     AliDebug(1,Form("TOF Window max dy: %f",dyMax));
 525     AliDebug(1,Form("TOF Window max dz: %f",dzMax));
 526     AliDebug(1,Form("TOF distance Cut: %f",dCut));
 527     AliDebug(1,Form("TOF Max Chi2: %f",maxChi2));
 528     AliDebug(1,Form("Time Walk Correction? : %d",timeWalkCorr));
 529   }
 530
 531   //Match ESD tracks to clusters in TOF
 532
 533   // Get the number of propagation steps
 534   Int_t nSteps=(Int_t)(detDepth/stepSize);
 535   AliDebug(1,Form(" Number of steps to be done %d",nSteps));
 536
 537   AliDebug(1,"++++++++++++++++++++++++++++++++++++++++++++++++++++++++");
 538
 539   //PH Arrays (moved outside of the loop)
 540   Float_t * trackPos[4];
 541   for (Int_t ii=0; ii<4; ii++) trackPos[ii] = new Float_t[nSteps];
 542   Int_t * clind = new Int_t[fN];
 543
 544   // Some init
 545   const Int_t kNclusterMax = 1000; // related to fN value
 546   TGeoHMatrix global[kNclusterMax];
 547
 548   //The matching loop
 549   for (Int_t iseed=0; iseed<fNseedsTOF; iseed++) {
 550
 551     fTOFtrackPoints->Delete();
 552
 553     for (Int_t ii=0; ii<kNclusterMax; ii++)
 554       global[ii] = 0x0;
 555     AliTOFtrack *track =(AliTOFtrack*)fTracks->UncheckedAt(iseed);
 556     AliESDtrack *t =(AliESDtrack*)fSeeds->At(track->GetSeedIndex());
 557     //if ( t->GetTOFsignal()>0. ) continue;
 558     if ( ((t->GetStatus()&AliESDtrack::kTOFout)!=0 ) && mLastStep < 2) continue;
 559     AliTOFtrack *trackTOFin = new AliTOFtrack(*track);
 560
 561     // Determine a window around the track
 562     Double_t x,par[5];
 563     trackTOFin->GetExternalParameters(x,par);
 564     Double_t cov[15];
 565     trackTOFin->GetExternalCovariance(cov);
 566
 567     if (cov[0]<0. || cov[2]<0.) {
 568       AliWarning(Form("Very strange track (%d)! At least one of its covariance matrix diagonal elements is negative!",iseed));
 569       //delete trackTOFin;
 570       //continue;
 571     }
 572
 573     Double_t dphi=
 574       scaleFact*
 575       ((5*TMath::Sqrt(TMath::Abs(cov[0])) + 0.5*dY + 2.5*TMath::Abs(par[2]))/sensRadius);
 576     Double_t dz=
 577        scaleFact*
 578        (5*TMath::Sqrt(TMath::Abs(cov[2])) + 0.5*dZ + 2.5*TMath::Abs(par[3]));
 579
 580     Double_t phi=TMath::ATan2(par[0],x) + trackTOFin->GetAlpha();
 581     if (phi<-TMath::Pi())phi+=2*TMath::Pi();
 582     if (phi>=TMath::Pi())phi-=2*TMath::Pi();
 583     Double_t z=par[1];
 584
 585     //upper limit on window's size.
 586     if (dz> dzMax) dz=dzMax;
 587     if (dphi*sensRadius> dyMax) dphi=dyMax/sensRadius;
 588
 589
 590     // find the clusters in the window of the track
 591     Int_t nc=0;
 592     for (Int_t k=FindClusterIndex(z-dz); k<fN; k++) {
 593
 594       if (nc>=kNclusterMax) {
 595         AliWarning("No more matchable clusters can be stored! Please, increase the corresponding vectors size.");
 596         break;
 597       }
 598
 599       AliTOFcluster *c=fClusters[k];
 600       if (c->GetZ() > z+dz) break;
 601       if (c->IsUsed() && mLastStep < 2) continue;
 602       if (!c->GetStatus()) {
 603         AliDebug(1,"Cluster in channel declared bad!");
 604         continue; // skip bad channels as declared in OCDB
 605       }
 606
 607       Double_t dph=TMath::Abs(c->GetPhi()-phi);
 608       if (dph>TMath::Pi()) dph-=2.*TMath::Pi();
 609       if (TMath::Abs(dph)>dphi) continue;
 610
 611       Double_t yc=(c->GetPhi() - trackTOFin->GetAlpha())*c->GetR();
 612       Double_t p[2]={yc, c->GetZ()};
 613       Double_t cov2[3]= {dY*dY/12., 0., dZ*dZ/12.};
 614       if (trackTOFin->AliExternalTrackParam::GetPredictedChi2(p,cov2) > maxChi2)continue;
 615
 616       clind[nc] = k;
 617       Char_t path[200];
 618       Int_t ind[5];
 619       ind[0]=c->GetDetInd(0);
 620       ind[1]=c->GetDetInd(1);
 621       ind[2]=c->GetDetInd(2);
 622       ind[3]=c->GetDetInd(3);
 623       ind[4]=c->GetDetInd(4);
 624       fGeom->GetVolumePath(ind,path);
 625       gGeoManager->cd(path);
 626       global[nc] = *gGeoManager->GetCurrentMatrix();
 627       nc++;
 628     }
 629
 630     if (nc == 0 ) {
 631       AliDebug(1,Form("No available clusters for the track number %d",iseed));
 632       fnunmatch++;
 633       delete trackTOFin;
 634       continue;
 635     }
 636
 637     AliDebug(1,Form(" Number of available TOF clusters for the track number %d: %d",iseed,nc));
 638
 639     //start fine propagation
 640
 641     Int_t nStepsDone = 0;
 642     for( Int_t istep=0; istep<nSteps; istep++){
 643
 644       // First of all, propagate the track...
 645       Float_t xs = AliTOFGeometry::RinTOF()+istep*stepSize;
 646       if (!(trackTOFin->PropagateTo(xs))) break;
 647
 648       //  ...and then, if necessary, rotate the track
 649       Double_t ymax = xs*TMath::Tan(0.5*AliTOFGeometry::GetAlpha());
 650       Double_t ysect = trackTOFin->GetY();
 651       if (ysect > ymax) {
 652         if (!(trackTOFin->Rotate(AliTOFGeometry::GetAlpha()))) break;
 653       } else if (ysect <-ymax) {
 654         if (!(trackTOFin->Rotate(-AliTOFGeometry::GetAlpha()))) break;
 655       }
 656
 657       nStepsDone++;
 658       AliDebug(3,Form(" current step %d (%d) - nStepsDone=%d",istep,nSteps,nStepsDone));
 659
 660       // store the running point (Globalrf) - fine propagation
 661
 662       Double_t r[3];
 663       trackTOFin->GetXYZ(r);
 664       trackPos[0][nStepsDone-1]= (Float_t) r[0];
 665       trackPos[1][nStepsDone-1]= (Float_t) r[1];
 666       trackPos[2][nStepsDone-1]= (Float_t) r[2];
 667       trackPos[3][nStepsDone-1]= trackTOFin->GetIntegratedLength();
 668     }
 669
 670
 671 #if 0
 672     /*****************/
 673     /**** OLD CODE ***/
 674     /*****************/
 675
 676     Int_t nfound = 0;
 677     Bool_t accept = kFALSE;
 678     Bool_t isInside = kFALSE;
 679     for (Int_t istep=0; istep<nStepsDone; istep++) {
 680
 681       Float_t ctrackPos[3];
 682       ctrackPos[0] = trackPos[0][istep];
 683       ctrackPos[1] = trackPos[1][istep];
 684       ctrackPos[2] = trackPos[2][istep];
 685
 686       //now see whether the track matches any of the TOF clusters
 687
 688       Float_t dist3d[3];
 689       accept = kFALSE;
 690       for (Int_t i=0; i<nc; i++) {
 691         isInside = fGeom->IsInsideThePad((TGeoHMatrix*)(&global[i]),ctrackPos,dist3d);
 692
 693         if ( mLastStep ) {
 694           Float_t yLoc = dist3d[1];
 695           Float_t rLoc = TMath::Sqrt(dist3d[0]*dist3d[0]+dist3d[2]*dist3d[2]);
 696           accept = (TMath::Abs(yLoc)<padDepth*0.5 && rLoc<dCut);
 697           AliDebug(3," I am in the case mLastStep==kTRUE ");
 698         }
 699         else {
 700           accept = isInside;
 701         }
 702         if (accept) {
 703
 704           fTOFtrackPoints->AddLast(new AliTOFtrackPoint(clind[i],
 705                                                         TMath::Sqrt(dist3d[0]*dist3d[0] + dist3d[1]*dist3d[1] + dist3d[2]*dist3d[2]),
 706                                                         dist3d[2],dist3d[0],dist3d[1],
 707                                                         AliTOFGeometry::RinTOF()+istep*stepSize,trackPos[3][istep]));
 708
 709           AliDebug(3,Form(" dist3dLoc[0] = %f, dist3dLoc[1] = %f, dist3dLoc[2] = %f ",dist3d[0],dist3d[1],dist3d[2]));
 710           nfound++;
 711           if(accept &&!mLastStep)break;
 712         }//end if accept
 713
 714       } //end for on the clusters
 715       if(accept &&!mLastStep)break;
 716     } //end for on the steps
 717
 718     /*****************/
 719     /**** OLD CODE ***/
 720     /*****************/
 721 #endif
 722
 723     if ( nStepsDone == 0 ) {
 724       AliDebug(1,Form(" No track points for the track number %d",iseed));
 725       fnunmatch++;
 726       delete trackTOFin;
 727       continue;
 728     }
 729
 730     AliDebug(2,Form(" Number of steps done for the track number %d: %d",iseed,nStepsDone));
 731
 732     /*****************/
 733     /**** NEW CODE ***/
 734     /*****************/
 735
 736     Int_t *isClusterMatchable = NULL;
 737     if(nc){
 738       isClusterMatchable = new Int_t[nc];
 739       for (Int_t i=0; i<nc; i++) isClusterMatchable[i] = kFALSE;
 740     }
 741
 742     Int_t nfound = 0;
 743     Bool_t accept = kFALSE;
 744     Bool_t isInside = kFALSE;
 745     for (Int_t istep=0; istep<nStepsDone; istep++) {
 746
 747       Bool_t gotInsideCluster = kFALSE;
 748       Int_t trackInsideCluster = -1;
 749
 750       Float_t ctrackPos[3];
 751       ctrackPos[0] = trackPos[0][istep];
 752       ctrackPos[1] = trackPos[1][istep];
 753       ctrackPos[2] = trackPos[2][istep];
 754
 755       //now see whether the track matches any of the TOF clusters
 756
 757       Float_t dist3d[3]={0.,0.,0.};
 758       accept = kFALSE;
 759       for (Int_t i=0; i<nc; i++) {
 760
 761         // ***** NEW *****
 762         /* check whether track was inside another cluster
 763          * and in case inhibit this cluster.
 764          * this will allow to only go on and add track points for
 765          * that cluster where the track got inside first */
 766         if (gotInsideCluster && trackInsideCluster != i) {
 767           AliDebug(3,Form(" A - istep=%d ~ %d %d ~ nfound=%d",istep,trackInsideCluster,i,nfound));
 768           continue;
 769         }
 770         AliDebug(3,Form(" B - istep=%d ~ %d %d ~ nfound=%d",istep,trackInsideCluster,i,nfound));
 771
 772         /* check whether track is inside this cluster */
 773         for (Int_t hh=0; hh<3; hh++) dist3d[hh]=0.;
 774         isInside = fGeom->IsInsideThePad((TGeoHMatrix*)(&global[i]),ctrackPos,dist3d);
 775
 776         // ***** NEW *****
 777         /* if track is inside this cluster set flags which will then
 778          * inhibit to add track points for the other clusters */
 779         if (isInside) {
 780           gotInsideCluster = kTRUE;
 781           trackInsideCluster = i;
 782         }
 783
 784         if ( mLastStep ) {
 785           Float_t yLoc = dist3d[1];
 786           Float_t rLoc = TMath::Sqrt(dist3d[0]*dist3d[0]+dist3d[2]*dist3d[2]);
 787           accept = (TMath::Abs(yLoc)<padDepth*0.5 && rLoc<dCut);
 788           AliDebug(3," I am in the case mLastStep==kTRUE ");
 789         }
 790
 791         //***** NEW *****
 792         /* add point everytime that:
 793          * - the track is inside the cluster
 794          * - the track got inside the cluster, even when it eventually exited the cluster
 795          * - the tracks is within dCut from the cluster
 796          */
 797         if (accept || isInside || gotInsideCluster) {
 798
 799           fTOFtrackPoints->AddLast(new AliTOFtrackPoint(clind[i],
 800                                                         TMath::Sqrt(dist3d[0]*dist3d[0] + dist3d[1]*dist3d[1] + dist3d[2]*dist3d[2]),
 801                                                         dist3d[2],dist3d[0],dist3d[1],
 802                                                         AliTOFGeometry::RinTOF()+istep*stepSize,trackPos[3][istep]));
 803
 804           AliDebug(2,Form(" dist3dLoc[0] = %f, dist3dLoc[1] = %f, dist3dLoc[2] = %f ",dist3d[0],dist3d[1],dist3d[2]));
 805           nfound++;
 806
 807           AliDebug(3,Form(" C - istep=%d ~ %d %d ~ nfound=%d",istep,trackInsideCluster,i,nfound));
 808
 809           // store the match in the ESD
 810           if (mLastStep==2 && !isClusterMatchable[i]) { // add TOF clusters to the track
 811             //
 812             isClusterMatchable[i] = kTRUE;
 813             //Tracking info
 814             Double_t mom=t->GetP();
 815             AliDebug(3,Form(" Momentum for track %d -> %f", iseed,mom));
 816             Double_t time[AliPID::kSPECIESC];
 817             // read from old structure (the one used by TPC in reco)
 818             for(Int_t isp=0;isp<AliPID::kSPECIESC;isp++){
 819               time[isp] = t->GetIntegratedTimesOld(isp); // in ps
 820               Double_t mass=AliPID::ParticleMass(isp);
 821               Double_t momz = mom*AliPID::ParticleCharge(isp);
 822               time[isp]+=(trackPos[3][istep]-trackPos[3][0])/kSpeedOfLight*TMath::Sqrt(momz*momz+mass*mass)/momz;
 823               //time[isp]+=(trackPos[3][istep]-trackPos[3][0])/kSpeedOfLight*TMath::Sqrt(mom*mom+mass*mass)/mom;
 824             }
 825             //
 826             AliESDTOFCluster* esdTOFCl = GetESDTOFCluster(clind[i]);
 827             if(!esdTOFCl->Update(t->GetID(),dist3d[0],dist3d[1],dist3d[2],trackPos[3][istep],time))//x,y,z -> tracking RF
 828               t->AddTOFcluster(esdTOFCl->GetESDID());
 829           }
 830
 831           // ***** NEW *****
 832           /* do not break loop in any case
 833            * if the track got inside a cluster all other clusters
 834            * are inhibited */
 835           //      if(accept &&!mLastStep)break;
 836
 837         }//end if accept
 838
 839       } //end for on the clusters
 840
 841       // ***** NEW *****
 842       /* do not break loop in any case
 843       * if the track got inside a cluster all other clusters
 844       * are inhibited but we want to go on adding track points */
 845       //      if(accept &&!mLastStep)break;
 846
 847     } //end for on the steps
 848     if(nc) delete[] isClusterMatchable;
 849
 850     if (nfound == 0 ) {
 851       AliDebug(1,Form("No track points for the track number %d",iseed));
 852       fnunmatch++;
 853       delete trackTOFin;
 854       continue;
 855     }
 856
 857     AliDebug(1,Form(" Number of track points for the track number %d: %d",iseed,nfound));
 858
 859     // now choose the cluster to be matched with the track.
 860
 861     Int_t idclus=-1;
 862     Float_t  recL = 0.;
 863     Float_t  xpos=0.;
 864     Float_t  mindist=1000.;
 865     Float_t  mindistZ=0.;
 866     Float_t  mindistY=0.;
 867     Float_t  mindistX=stepSize;
 868     for (Int_t iclus= 0; iclus<nfound;iclus++) {
 869       AliTOFtrackPoint *matchableTOFcluster = (AliTOFtrackPoint*)fTOFtrackPoints->At(iclus);
 870       //if ( matchableTOFcluster->Distance()<mindist ) {
 871       if ( TMath::Abs(matchableTOFcluster->DistanceX())<TMath::Abs(mindistX) &&
 872            TMath::Abs(matchableTOFcluster->DistanceX())<=stepSize ) {
 873         mindist = matchableTOFcluster->Distance();
 874         mindistZ = matchableTOFcluster->DistanceZ(); // Z distance in the
 875                                                      // RF of the hit pad
 876                                                      // closest to the
 877                                                      // reconstructed
 878                                                      // track
 879         mindistY = matchableTOFcluster->DistanceY(); // Y distance in the
 880                                                      // RF of the hit pad
 881                                                      // closest to the
 882                                                      // reconstructed
 883                                                      // track
 884         mindistX = matchableTOFcluster->DistanceX(); // X distance in the
 885                                                      // RF of the hit pad
 886                                                      // closest to the
 887                                                      // reconstructed
 888                                                      // track
 889         xpos = matchableTOFcluster->PropRadius();
 890         idclus = matchableTOFcluster->Index();
 891         recL = matchableTOFcluster->Length();// + corrLen*0.5;
 892
 893         AliDebug(2,Form(" %d(%d) --- %f (%f, %f, %f), step=%f -- idclus=%d --- seed=%d, trackId=%d, trackLab=%d", iclus,nfound,
 894                         mindist,mindistX,mindistY,mindistZ,stepSize,idclus,iseed,track->GetSeedIndex(),track->GetLabel()));
 895
 896       }
 897     } // loop on found TOF track points
 898
 899     if (TMath::Abs(mindistX)>stepSize && idclus!=-1) {
 900       AliInfo(Form(" %d - not matched --- but idclus=%d, trackId=%d, trackLab=%d",iseed,
 901                    idclus,track->GetSeedIndex(),track->GetLabel()));
 902       idclus=-1;
 903     }
 904
 905     if (idclus==-1) {
 906       AliDebug(1,Form("Reconstructed track %d doesn't match any TOF cluster", iseed));
 907       fnunmatch++;
 908       delete trackTOFin;
 909       continue;
 910     }
 911
 912     AliDebug(1,Form(" %d - matched",iseed));
 913
 914     fnmatch++;
 915
 916     AliTOFcluster *c=fClusters[idclus];
 917
 918     AliDebug(3, Form("%7d     %7d     %10d     %10d  %10d  %10d      %7d",
 919                      iseed,
 920                      fnmatch-1,
 921                      TMath::Abs(trackTOFin->GetLabel()),
 922                      c->GetLabel(0), c->GetLabel(1), c->GetLabel(2),
 923                      idclus)); // AdC
 924
 925     c->Use();
 926
 927     // Track length correction for matching Step 2
 928     /*
 929     if (mLastStep) {
 930       Float_t rc = TMath::Sqrt(c->GetR()*c->GetR() + c->GetZ()*c->GetZ());
 931       Float_t rt = TMath::Sqrt(trackPos[0][70]*trackPos[0][70]
 932                                +trackPos[1][70]*trackPos[1][70]
 933                                +trackPos[2][70]*trackPos[2][70]);
 934       Float_t dlt=rc-rt;
 935       recL=trackPos[3][70]+dlt;
 936     }
 937     */
 938     if (
 939         (c->GetLabel(0)==TMath::Abs(trackTOFin->GetLabel()))
 940         ||
 941         (c->GetLabel(1)==TMath::Abs(trackTOFin->GetLabel()))
 942         ||
 943         (c->GetLabel(2)==TMath::Abs(trackTOFin->GetLabel()))
 944         ) {
 945       fngoodmatch++;
 946
 947        AliDebug(2,Form(" track label good %5d",trackTOFin->GetLabel()));
 948
 949     }
 950     else {
 951       fnbadmatch++;
 952
 953       AliDebug(2,Form(" track label  bad %5d",trackTOFin->GetLabel()));
 954
 955     }
 956
 957     delete trackTOFin;
 958
 959     //  Store quantities to be used in the TOF Calibration
 960     Float_t tToT=AliTOFGeometry::ToTBinWidth()*c->GetToT()*1E-3; // in ns
 961 //     t->SetTOFsignalToT(tToT);
 962     Float_t rawTime=AliTOFGeometry::TdcBinWidth()*c->GetTDCRAW()+kTimeOffset; // RAW time,in ps
 963 //     t->SetTOFsignalRaw(rawTime);
 964 //     t->SetTOFsignalDz(mindistZ);
 965 //     t->SetTOFsignalDx(mindistY);
 966 //     t->SetTOFDeltaBC(c->GetDeltaBC());
 967 //     t->SetTOFL0L1(c->GetL0L1Latency());
 968
 969     Float_t info[10] = {mindist,mindistY,mindistZ,
 970                         0.,0.,0.,0.,0.,0.,0.};
 971     t->SetTOFInfo(info);
 972     AliDebug(3,Form(" distance=%f; residual in the pad reference frame: dX=%f, dZ=%f", info[0],info[1],info[2]));
 973
 974
 975     Int_t ind[5];
 976     ind[0]=c->GetDetInd(0);
 977     ind[1]=c->GetDetInd(1);
 978     ind[2]=c->GetDetInd(2);
 979     ind[3]=c->GetDetInd(3);
 980     ind[4]=c->GetDetInd(4);
 981     Int_t calindex = AliTOFGeometry::GetIndex(ind);
 982 //     t->SetTOFCalChannel(calindex);
 983
 984     // keep track of the track labels in the matched cluster
 985     Int_t tlab[3];
 986     tlab[0]=c->GetLabel(0);
 987     tlab[1]=c->GetLabel(1);
 988     tlab[2]=c->GetLabel(2);
 989     AliDebug(3,Form(" tdc time of the matched track %6d = ",c->GetTDC()));
 990     Double_t tof=AliTOFGeometry::TdcBinWidth()*c->GetTDC()+kTimeOffset; // in ps
 991     AliDebug(3,Form(" tof time of the matched track: %f = ",tof));
 992     Double_t tofcorr=tof;
 993     if(timeWalkCorr)tofcorr=CorrectTimeWalk(mindistZ,tof);
 994     AliDebug(3,Form(" tof time of the matched track, after TW corr: %f = ",tofcorr));
 995     //Set TOF time signal and pointer to the matched cluster
 996 //     t->SetTOFsignal(tofcorr);
 997     t->SetTOFcluster(idclus); // pointing to the recPoints tree
 998
 999     AliDebug(3,Form(" Setting TOF raw time: %f, z distance: %f  corrected time: %f ",rawTime,mindistZ,tofcorr));
1000
1001     //Tracking info
1002     Double_t time[AliPID::kSPECIESC];
1003     // read from old structure (the one used by TPC in reco)
1004     for(Int_t isp=0;isp<AliPID::kSPECIESC;isp++){
1005       time[isp] = t->GetIntegratedTimesOld(isp); // in ps
1006     }
1007     Double_t mom=t->GetP();
1008     AliDebug(3,Form(" Momentum for track %d -> %f", iseed,mom));
1009     for (Int_t j=0;j<AliPID::kSPECIESC;j++) {
1010       Double_t mass=AliPID::ParticleMass(j);
1011       Double_t momz = mom*AliPID::ParticleCharge(j);
1012       time[j]+=(recL-trackPos[3][0])/kSpeedOfLight*TMath::Sqrt(momz*momz+mass*mass)/momz;
1013       //time[j]+=(recL-trackPos[3][0])/kSpeedOfLight*TMath::Sqrt(mom*mom+mass*mass)/mom;
1014     }
1015
1016     AliTOFtrack *trackTOFout = new AliTOFtrack(*t);
1017     if (!(trackTOFout->PropagateTo(xpos))) {
1018       delete trackTOFout;
1019       continue;
1020     }
1021
1022     // If necessary, rotate the track
1023     Double_t yATxposMax=xpos*TMath::Tan(0.5*AliTOFGeometry::GetAlpha());
1024     Double_t yATxpos=trackTOFout->GetY();
1025     if (yATxpos > yATxposMax) {
1026       if (!(trackTOFout->Rotate(AliTOFGeometry::GetAlpha()))) {
1027         delete trackTOFout;
1028         continue;
1029       }
1030     } else if (yATxpos < -yATxposMax) {
1031       if (!(trackTOFout->Rotate(-AliTOFGeometry::GetAlpha()))) {
1032         delete trackTOFout;
1033         continue;
1034       }
1035     }
1036
1037     // Fill the track residual histograms and update track only if in the first two step (0 and 1)
1038     if(mLastStep < 2){
1039       FillResiduals(trackTOFout,c,kFALSE);
1040
1041       t->UpdateTrackParams(trackTOFout,AliESDtrack::kTOFout);
1042
1043 // don't update old structure with TOF info
1044 //       t->SetIntegratedLength(recL);
1045 //       t->SetIntegratedTimes(time);
1046 //       t->SetTOFLabel(tlab);
1047
1048       // add tof cluster to the track also for step 2
1049       AliESDTOFCluster* esdTOFCl = GetESDTOFCluster(idclus);
1050       esdTOFCl->Update(t->GetID(),mindistY,mindist,mindistZ,recL,time);
1051       t->AddTOFcluster(esdTOFCl->GetESDID());
1052       /* RS?
1053       if(idclus < 20000){
1054         fClusterESD[idclus]->Update(t->GetID(),mindistY,mindist,mindistZ,recL,time);//x,y,z -> tracking RF
1055
1056         t->AddTOFcluster(idclus);
1057       }
1058       else{
1059         AliInfo("Too many TOF clusters matched with tracks (> 20000)");
1060       }
1061       */
1062     }
1063     // Fill Reco-QA histos for Reconstruction
1064     fHRecNClus->Fill(nc);
1065     fHRecDist->Fill(mindist);
1066     if (cov[0]>=0.)
1067       fHRecSigYVsP->Fill(mom,TMath::Sqrt(cov[0]));
1068     else
1069       fHRecSigYVsP->Fill(mom,-TMath::Sqrt(-cov[0]));
1070     if (cov[2]>=0.)
1071       fHRecSigZVsP->Fill(mom,TMath::Sqrt(cov[2]));
1072     else
1073       fHRecSigZVsP->Fill(mom,-TMath::Sqrt(-cov[2]));
1074     fHRecSigYVsPWin->Fill(mom,dphi*sensRadius);
1075     fHRecSigZVsPWin->Fill(mom,dz);
1076
1077     // Fill Tree for on-the-fly offline Calibration
1078
1079     if ( !((t->GetStatus() & AliESDtrack::kTIME)==0 ) ) {
1080       fIch=calindex;
1081       fToT=tToT;
1082       fTime=rawTime;
1083       fExpTimePi=time[2];
1084       fExpTimeKa=time[3];
1085       fExpTimePr=time[4];
1086       fCalTree->Fill();
1087     }
1088     delete trackTOFout;
1089   }
1090
1091
1092   for (Int_t ii=0; ii<4; ii++) delete [] trackPos[ii];
1093   delete [] clind;
1094
1095 }
1096 //_________________________________________________________________________
1097 Int_t AliTOFtracker::LoadClusters(TTree *cTree) {
1098   //--------------------------------------------------------------------
1099   //This function loads the TOF clusters
1100   //--------------------------------------------------------------------
1101
1102   Int_t npadX = AliTOFGeometry::NpadX();
1103   Int_t npadZ = AliTOFGeometry::NpadZ();
1104   Int_t nStripA = AliTOFGeometry::NStripA();
1105   Int_t nStripB = AliTOFGeometry::NStripB();
1106   Int_t nStripC = AliTOFGeometry::NStripC();
1107
1108   TBranch *branch=cTree->GetBranch("TOF");
1109   if (!branch) {
1110     AliError("can't get the branch with the TOF clusters !");
1111     return 1;
1112   }
1113
1114   static TClonesArray dummy("AliTOFcluster",10000);
1115   dummy.Clear();
1116   TClonesArray *clusters=&dummy;
1117   branch->SetAddress(&clusters);
1118
1119   cTree->GetEvent(0);
1120   Int_t nc=clusters->GetEntriesFast();
1121   fHDigNClus->Fill(nc);
1122
1123   AliInfo(Form("Number of clusters: %d",nc));
1124
1125   fN = 0;
1126   fNTOFmatched = 0;
1127
1128   for (Int_t i=0; i<nc; i++) {
1129     AliTOFcluster *c=(AliTOFcluster*)clusters->UncheckedAt(i);
1130 //PH    fClusters[i]=new AliTOFcluster(*c); fN++;
1131     fClusters[i]=c; fN++;
1132     c->SetESDID(-1);
1133   // Fill Digits QA histos
1134
1135     Int_t isector = c->GetDetInd(0);
1136     Int_t iplate = c->GetDetInd(1);
1137     Int_t istrip = c->GetDetInd(2);
1138     Int_t ipadX = c->GetDetInd(4);
1139     Int_t ipadZ = c->GetDetInd(3);
1140
1141     Float_t time =(AliTOFGeometry::TdcBinWidth()*c->GetTDC())*1E-3; // in ns
1142     Float_t tot = (AliTOFGeometry::TdcBinWidth()*c->GetToT())*1E-3;//in ns
1143
1144     /* RS?
1145     Int_t ind[5];
1146     ind[0]=isector;
1147     ind[1]=iplate;
1148     ind[2]=istrip;
1149     ind[3]=ipadZ;
1150     ind[4]=ipadX;
1151     Int_t calindex = AliTOFGeometry::GetIndex(ind);
1152     Int_t tofLabels[3]={c->GetLabel(0),c->GetLabel(1),c->GetLabel(2)};
1153     */
1154     Int_t stripOffset = 0;
1155     switch (iplate) {
1156     case 0:
1157       stripOffset = 0;
1158       break;
1159     case 1:
1160       stripOffset = nStripC;
1161       break;
1162     case 2:
1163       stripOffset = nStripC+nStripB;
1164       break;
1165     case 3:
1166       stripOffset = nStripC+nStripB+nStripA;
1167       break;
1168     case 4:
1169       stripOffset = nStripC+nStripB+nStripA+nStripB;
1170       break;
1171     default:
1172       AliError(Form("Wrong plate number in TOF (%d) !",iplate));
1173       break;
1174     };
1175     Int_t zindex=npadZ*(istrip+stripOffset)+(ipadZ+1);
1176     Int_t phiindex=npadX*isector+ipadX+1;
1177     fHDigClusMap->Fill(zindex,phiindex);
1178     fHDigClusTime->Fill(time);
1179     fHDigClusToT->Fill(tot);
1180
1181     fNTOFmatched++; // RS: Actually number of clusters
1182     /* RS?
1183     if(fNTOFmatched < 20000){
1184       fHit[fNTOFmatched] = new AliESDTOFHit(AliTOFGeometry::TdcBinWidth()*c->GetTDC(),
1185                                             AliTOFGeometry::TdcBinWidth()*c->GetTDCRAW(),
1186                                             AliTOFGeometry::ToTBinWidth()*c->GetToT()*1E-3,
1187                                             calindex,tofLabels,c->GetL0L1Latency(),
1188                                             c->GetDeltaBC(),i,c->GetZ(),c->GetR(),c->GetPhi());
1189       fNTOFmatched++;
1190     }
1191     */
1192   }
1193
1194   return 0;
1195 }
1196 //_________________________________________________________________________
1197 void AliTOFtracker::UnloadClusters() {
1198   //--------------------------------------------------------------------
1199   //This function unloads TOF clusters
1200   //--------------------------------------------------------------------
1201   for (Int_t i=0; i<fN; i++) {
1202 //PH    delete fClusters[i];
1203     fClusters[i] = 0x0;
1204   }
1205   /* RS
1206   for(Int_t i=0; i< 20000;i++){
1207     if(fClusterESD[i]){
1208       delete fClusterESD[i];
1209       fClusterESD[i] = NULL;
1210     }
1211     if(fHit[i]){
1212       delete fHit[i];
1213       fHit[i] = NULL;
1214     }
1215   }
1216   */
1217   fN=0;
1218   fNTOFmatched = 0;
1219 }
1220
1221 //_________________________________________________________________________
1222 Int_t AliTOFtracker::FindClusterIndex(Double_t z) const {
1223   //--------------------------------------------------------------------
1224   // This function returns the index of the nearest cluster
1225   //--------------------------------------------------------------------
1226   if (fN==0) return 0;
1227   if (z <= fClusters[0]->GetZ()) return 0;
1228   if (z > fClusters[fN-1]->GetZ()) return fN;
1229   Int_t b=0, e=fN-1, m=(b+e)/2;
1230   for (; b<e; m=(b+e)/2) {
1231     if (z > fClusters[m]->GetZ()) b=m+1;
1232     else e=m;
1233   }
1234   return m;
1235 }
1236
1237 //_________________________________________________________________________
1238 Bool_t AliTOFtracker::GetTrackPoint(Int_t index, AliTrackPoint& p) const
1239 {
1240   // Get track space point with index i
1241   // Coordinates are in the global system
1242   AliTOFcluster *cl = fClusters[index];
1243   Float_t xyz[3];
1244   xyz[0] = cl->GetR()*TMath::Cos(cl->GetPhi());
1245   xyz[1] = cl->GetR()*TMath::Sin(cl->GetPhi());
1246   xyz[2] = cl->GetZ();
1247   Float_t phiangle = (Int_t(cl->GetPhi()*TMath::RadToDeg()/20.)+0.5)*20.*TMath::DegToRad();
1248   Float_t sinphi = TMath::Sin(phiangle), cosphi = TMath::Cos(phiangle);
1249   Float_t tiltangle = AliTOFGeometry::GetAngles(cl->GetDetInd(1),cl->GetDetInd(2))*TMath::DegToRad();
1250   Float_t sinth = TMath::Sin(tiltangle), costh = TMath::Cos(tiltangle);
1251   Float_t sigmay2 = AliTOFGeometry::XPad()*AliTOFGeometry::XPad()/12.;
1252   Float_t sigmaz2 = AliTOFGeometry::ZPad()*AliTOFGeometry::ZPad()/12.;
1253   Float_t cov[6];
1254   cov[0] = sinphi*sinphi*sigmay2 + cosphi*cosphi*sinth*sinth*sigmaz2;
1255   cov[1] = -sinphi*cosphi*sigmay2 + sinphi*cosphi*sinth*sinth*sigmaz2;
1256   cov[2] = -cosphi*sinth*costh*sigmaz2;
1257   cov[3] = cosphi*cosphi*sigmay2 + sinphi*sinphi*sinth*sinth*sigmaz2;
1258   cov[4] = -sinphi*sinth*costh*sigmaz2;
1259   cov[5] = costh*costh*sigmaz2;
1260   p.SetXYZ(xyz[0],xyz[1],xyz[2],cov);
1261
1262   // Detector numbering scheme
1263   Int_t nSector = AliTOFGeometry::NSectors();
1264   Int_t nPlate  = AliTOFGeometry::NPlates();
1265   Int_t nStripA = AliTOFGeometry::NStripA();
1266   Int_t nStripB = AliTOFGeometry::NStripB();
1267   Int_t nStripC = AliTOFGeometry::NStripC();
1268
1269   Int_t isector = cl->GetDetInd(0);
1270   if (isector >= nSector)
1271     AliError(Form("Wrong sector number in TOF (%d) !",isector));
1272   Int_t iplate = cl->GetDetInd(1);
1273   if (iplate >= nPlate)
1274     AliError(Form("Wrong plate number in TOF (%d) !",iplate));
1275   Int_t istrip = cl->GetDetInd(2);
1276
1277   Int_t stripOffset = 0;
1278   switch (iplate) {
1279   case 0:
1280     stripOffset = 0;
1281     break;
1282   case 1:
1283     stripOffset = nStripC;
1284     break;
1285   case 2:
1286     stripOffset = nStripC+nStripB;
1287     break;
1288   case 3:
1289     stripOffset = nStripC+nStripB+nStripA;
1290     break;
1291   case 4:
1292     stripOffset = nStripC+nStripB+nStripA+nStripB;
1293     break;
1294   default:
1295     AliError(Form("Wrong plate number in TOF (%d) !",iplate));
1296     break;
1297   };
1298
1299   Int_t idet = (2*(nStripC+nStripB)+nStripA)*isector +
1300                stripOffset +
1301                istrip;
1302   UShort_t volid = AliGeomManager::LayerToVolUID(AliGeomManager::kTOF,idet);
1303   p.SetVolumeID((UShort_t)volid);
1304   return kTRUE;
1305 }
1306 //_________________________________________________________________________
1307 void AliTOFtracker::InitCheckHists() {
1308
1309   //Init histos for Digits/Reco QA and Calibration
1310
1311
1312   TDirectory *dir = gDirectory;
1313   TFile *logFileTOF = 0;
1314
1315   TSeqCollection *list = gROOT->GetListOfFiles();
1316   int n = list->GetEntries();
1317   Bool_t isThere=kFALSE;
1318   for(int i=0; i<n; i++) {
1319     logFileTOF = (TFile*)list->At(i);
1320     if (strstr(logFileTOF->GetName(), "TOFQA.root")){
1321       isThere=kTRUE;
1322       break;
1323     }
1324   }
1325
1326   if(!isThere)logFileTOF = new TFile( "TOFQA.root","RECREATE");
1327   logFileTOF->cd();
1328
1329   fCalTree = new TTree("CalTree", "Tree for TOF calibration");
1330   fCalTree->Branch("TOFchannelindex",&fIch,"iTOFch/I");
1331   fCalTree->Branch("ToT",&fToT,"TOFToT/F");
1332   fCalTree->Branch("TOFtime",&fTime,"TOFtime/F");
1333   fCalTree->Branch("PionExpTime",&fExpTimePi,"PiExpTime/F");
1334   fCalTree->Branch("KaonExpTime",&fExpTimeKa,"KaExpTime/F");
1335   fCalTree->Branch("ProtonExpTime",&fExpTimePr,"PrExpTime/F");
1336
1337   //Digits "QA"
1338   fHDigClusMap = new TH2F("TOFDig_ClusMap", "",182,0.5,182.5,864, 0.5,864.5);
1339   fHDigNClus = new TH1F("TOFDig_NClus", "",200,0.5,200.5);
1340   fHDigClusTime = new TH1F("TOFDig_ClusTime", "",2000,0.,200.);
1341   fHDigClusToT = new TH1F("TOFDig_ClusToT", "",500,0.,100);
1342
1343   //Reco "QA"
1344   fHRecNClus =new TH1F("TOFRec_NClusW", "",50,0.5,50.5);
1345   fHRecDist=new TH1F("TOFRec_Dist", "",50,0.5,10.5);
1346   fHRecSigYVsP=new TH2F("TOFDig_SigYVsP", "",40,0.,4.,100, 0.,5.);
1347   fHRecSigZVsP=new TH2F("TOFDig_SigZVsP", "",40,0.,4.,100, 0.,5.);
1348   fHRecSigYVsPWin=new TH2F("TOFDig_SigYVsPWin", "",40,0.,4.,100, 0.,50.);
1349   fHRecSigZVsPWin=new TH2F("TOFDig_SigZVsPWin", "",40,0.,4.,100, 0.,50.);
1350
1351   dir->cd();
1352
1353 }
1354
1355 //_________________________________________________________________________
1356 void AliTOFtracker::SaveCheckHists() {
1357
1358   //write histos for Digits/Reco QA and Calibration
1359
1360   TDirectory *dir = gDirectory;
1361   TFile *logFileTOF = 0;
1362
1363   TSeqCollection *list = gROOT->GetListOfFiles();
1364   int n = list->GetEntries();
1365   Bool_t isThere=kFALSE;
1366   for(int i=0; i<n; i++) {
1367     logFileTOF = (TFile*)list->At(i);
1368     if (strstr(logFileTOF->GetName(), "TOFQA.root")){
1369       isThere=kTRUE;
1370       break;
1371     }
1372   }
1373
1374   if(!isThere) {
1375           AliError(Form("File TOFQA.root not found!! not wring histograms...."));
1376           return;
1377   }
1378   logFileTOF->cd();
1379   fHDigClusMap->Write(fHDigClusMap->GetName(), TObject::kOverwrite);
1380   fHDigNClus->Write(fHDigNClus->GetName(), TObject::kOverwrite);
1381   fHDigClusTime->Write(fHDigClusTime->GetName(), TObject::kOverwrite);
1382   fHDigClusToT->Write(fHDigClusToT->GetName(), TObject::kOverwrite);
1383   fHRecNClus->Write(fHRecNClus->GetName(), TObject::kOverwrite);
1384   fHRecDist->Write(fHRecDist->GetName(), TObject::kOverwrite);
1385   fHRecSigYVsP->Write(fHRecSigYVsP->GetName(), TObject::kOverwrite);
1386   fHRecSigZVsP->Write(fHRecSigZVsP->GetName(), TObject::kOverwrite);
1387   fHRecSigYVsPWin->Write(fHRecSigYVsPWin->GetName(), TObject::kOverwrite);
1388   fHRecSigZVsPWin->Write(fHRecSigZVsPWin->GetName(), TObject::kOverwrite);
1389   fCalTree->Write(fCalTree->GetName(),TObject::kOverwrite);
1390   logFileTOF->Flush();
1391
1392   dir->cd();
1393   }
1394 //_________________________________________________________________________
1395 Float_t AliTOFtracker::CorrectTimeWalk( Float_t dist, Float_t tof) const {
1396
1397   //dummy, for the moment
1398   Float_t tofcorr=0.;
1399   if(dist<AliTOFGeometry::ZPad()*0.5){
1400     tofcorr=tof;
1401     //place here the actual correction
1402   }else{
1403     tofcorr=tof;
1404   }
1405   return tofcorr;
1406 }
1407 //_________________________________________________________________________
1408
1409 void AliTOFtracker::FillClusterArray(TObjArray* arr) const
1410 {
1411   //
1412   // Returns the TOF cluster array
1413   //
1414
1415   if (fN==0)
1416     arr = 0x0;
1417   else
1418     for (Int_t i=0; i<fN; ++i) arr->Add(fClusters[i]);
1419
1420 }
1421
1422 //_________________________________________________________________________
1423 AliESDTOFCluster* AliTOFtracker::GetESDTOFCluster(int clID)
1424 {
1425   // get ESDTOFcluster corresponding to fClusters[clID]. If the original cluster
1426   // was not stored yet in the ESD, first do this
1427   AliTOFcluster *c = fClusters[clID];
1428   AliESDTOFCluster *clESD = 0;
1429   int esdID = c->GetESDID(); // was this cluster already stored in the ESD clusters?
1430   TClonesArray* esdTOFClArr = fESDEv->GetESDTOFClusters();
1431   if (esdID<0) { // cluster was not stored yet, do this
1432     esdID = esdTOFClArr->GetEntriesFast();
1433     c->SetESDID(esdID);
1434     // first store the hits of the cluster
1435     TClonesArray* esdTOFHitArr = fESDEv->GetESDTOFHits();
1436     int nh = esdTOFHitArr->GetEntriesFast();
1437     Int_t tofLabels[3]={c->GetLabel(0),c->GetLabel(1),c->GetLabel(2)};
1438     Int_t ind[5] = {c->GetDetInd(0), c->GetDetInd(1), c->GetDetInd(2), c->GetDetInd(3), c->GetDetInd(4) };
1439     Int_t calindex = AliTOFGeometry::GetIndex(ind);
1440     /*AliESDTOFHit* esdHit = */
1441     new ( (*esdTOFHitArr)[nh] )
1442       AliESDTOFHit( AliTOFGeometry::TdcBinWidth()*c->GetTDC(),
1443                     AliTOFGeometry::TdcBinWidth()*c->GetTDCRAW(),
1444                     AliTOFGeometry::ToTBinWidth()*c->GetToT()*1E-3,
1445                     calindex,tofLabels,c->GetL0L1Latency(),
1446                     c->GetDeltaBC(),esdID,c->GetZ(),c->GetR(),c->GetPhi());
1447     //
1448     clESD =  new( (*esdTOFClArr)[esdID] ) AliESDTOFCluster( clID );
1449     clESD->SetEvent(fESDEv);
1450     clESD->SetStatus( c->GetStatus() );
1451     clESD->SetESDID(esdID);
1452     //
1453     // register hits in the cluster
1454     clESD->AddESDTOFHitIndex(nh);
1455   }
1456   else clESD =  (AliESDTOFCluster*)esdTOFClArr->At(esdID); // cluster is aready stored in the ESD
1457   //
1458   return clESD;
1459   //
1460 }