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   fTracks(new TClonesArray("AliTOFtrack")),
  74   fSeeds(new TObjArray(100)),
  75   fTOFtrackPoints(new TObjArray(10)),
  76   fHDigClusMap(0x0),
  77   fHDigNClus(0x0),
  78   fHDigClusTime(0x0),
  79   fHDigClusToT(0x0),
  80   fHRecNClus(0x0),
  81   fHRecDist(0x0),
  82   fHRecSigYVsP(0x0),
  83   fHRecSigZVsP(0x0),
  84   fHRecSigYVsPWin(0x0),
  85   fHRecSigZVsPWin(0x0),
  86   fCalTree(0x0),
  87   fIch(-1),
  88   fToT(-1.),
  89   fTime(-1.),
  90   fExpTimePi(-1.),
  91   fExpTimeKa(-1.),
  92   fExpTimePr(-1.)
  93  {
  94   //AliTOFtracker main Ctor
  95
  96    for (Int_t ii=0; ii<kMaxCluster; ii++) fClusters[ii]=0x0;
  97
  98    // Gettimg the geometry
  99    fGeom= new AliTOFGeometry();
 100
 101    InitCheckHists();
 102
 103 }
 104 //_____________________________________________________________________________
 105 AliTOFtracker::~AliTOFtracker() {
 106   //
 107   // Dtor
 108   //
 109
 110   SaveCheckHists();
 111
 112   if(!(AliCDBManager::Instance()->GetCacheFlag())){
 113     delete fkRecoParam;
 114   }
 115   delete fGeom;
 116   delete fHDigClusMap;
 117   delete fHDigNClus;
 118   delete fHDigClusTime;
 119   delete fHDigClusToT;
 120   delete fHRecNClus;
 121   delete fHRecDist;
 122   delete fHRecSigYVsP;
 123   delete fHRecSigZVsP;
 124   delete fHRecSigYVsPWin;
 125   delete fHRecSigZVsPWin;
 126   delete fCalTree;
 127   if (fTracks){
 128     fTracks->Delete();
 129     delete fTracks;
 130     fTracks=0x0;
 131   }
 132   if (fSeeds){
 133     fSeeds->Delete();
 134     delete fSeeds;
 135     fSeeds=0x0;
 136   }
 137   if (fTOFtrackPoints){
 138     fTOFtrackPoints->Delete();
 139     delete fTOFtrackPoints;
 140     fTOFtrackPoints=0x0;
 141   }
 142
 143   for (Int_t ii=0; ii<kMaxCluster; ii++)
 144     if (fClusters[ii]) fClusters[ii]->Delete();
 145
 146 }
 147 //_____________________________________________________________________________
 148 void AliTOFtracker::GetPidSettings(AliESDpid *esdPID) {
 149   //
 150   // Sets TOF resolution from RecoParams
 151   //
 152   if (fkRecoParam)
 153     esdPID->GetTOFResponse().SetTimeResolution(fkRecoParam->GetTimeResolution());
 154   else
 155     AliWarning("fkRecoParam not yet set; cannot set PID settings");
 156 }
 157 //_____________________________________________________________________________
 158 Int_t AliTOFtracker::PropagateBack(AliESDEvent * const event) {
 159   //
 160   // Gets seeds from ESD event and Match with TOF Clusters
 161   //
 162
 163   // initialize RecoParam for current event
 164   AliDebug(1,"Initializing params for TOF");
 165
 166   fkRecoParam = AliTOFReconstructor::GetRecoParam();  // instantiate reco param from STEER...
 167
 168   if (fkRecoParam == 0x0) {
 169     AliFatal("No Reco Param found for TOF!!!");
 170   }
 171   //fkRecoParam->Dump();
 172   //if(fkRecoParam->GetApplyPbPbCuts())fkRecoParam=fkRecoParam->GetPbPbparam();
 173   //fkRecoParam->PrintParameters();
 174
 175   //Initialise some counters
 176
 177   fNseeds=0;
 178   fNseedsTOF=0;
 179   fngoodmatch=0;
 180   fnbadmatch=0;
 181   fnunmatch=0;
 182   fnmatch=0;
 183
 184   Int_t ntrk=event->GetNumberOfTracks();
 185   fNseeds = ntrk;
 186
 187
 188   //Load ESD tracks into a local Array of ESD Seeds
 189   for (Int_t i=0; i<fNseeds; i++)
 190     fSeeds->AddLast(event->GetTrack(i));
 191
 192   //Prepare ESD tracks candidates for TOF Matching
 193   CollectESD();
 194
 195   //First Step with Strict Matching Criterion
 196   MatchTracks(kFALSE);
 197
 198   //Second Step with Looser Matching Criterion
 199   MatchTracks(kTRUE);
 200
 201   AliInfo(Form("Number of matched tracks = %d (good = %d, bad = %d)",fnmatch,fngoodmatch,fnbadmatch));
 202
 203   //Update the matched ESD tracks
 204
 205   for (Int_t i=0; i<ntrk; i++) {
 206     AliESDtrack *t=event->GetTrack(i);
 207     AliESDtrack *seed =(AliESDtrack*)fSeeds->At(i);
 208
 209     if ( (seed->GetStatus()&AliESDtrack::kTOFin)!=0 ) {
 210       t->SetStatus(AliESDtrack::kTOFin);
 211       //if(seed->GetTOFsignal()>0){
 212       if ( (seed->GetStatus()&AliESDtrack::kTOFout)!=0 ) {
 213         t->SetStatus(AliESDtrack::kTOFout);
 214         t->SetTOFsignal(seed->GetTOFsignal());
 215         t->SetTOFcluster(seed->GetTOFcluster());
 216         t->SetTOFsignalToT(seed->GetTOFsignalToT());
 217         t->SetTOFsignalRaw(seed->GetTOFsignalRaw());
 218         t->SetTOFsignalDz(seed->GetTOFsignalDz());
 219         t->SetTOFsignalDx(seed->GetTOFsignalDx());
 220         t->SetTOFDeltaBC(seed->GetTOFDeltaBC());
 221         t->SetTOFL0L1(seed->GetTOFL0L1());
 222         t->SetTOFCalChannel(seed->GetTOFCalChannel());
 223         Int_t tlab[3]; seed->GetTOFLabel(tlab);
 224         t->SetTOFLabel(tlab);
 225
 226         Double_t alphaA = (Double_t)t->GetAlpha();
 227         Double_t xA = (Double_t)t->GetX();
 228         Double_t yA = (Double_t)t->GetY();
 229         Double_t zA = (Double_t)t->GetZ();
 230         Double_t p1A = (Double_t)t->GetSnp();
 231         Double_t p2A = (Double_t)t->GetTgl();
 232         Double_t p3A = (Double_t)t->GetSigned1Pt();
 233         const Double_t *covA = (Double_t*)t->GetCovariance();
 234
 235         // Make attention, please:
 236         //      AliESDtrack::fTOFInfo array does not be stored in the AliESDs.root file
 237         //      it is there only for a check during the reconstruction step.
 238         Float_t info[10]; seed->GetTOFInfo(info);
 239         t->SetTOFInfo(info);
 240         AliDebug(3,Form(" distance=%f; residual in the pad reference frame: dX=%f, dZ=%f", info[0],info[1],info[2]));
 241
 242         // Check done:
 243         //       by calling the AliESDtrack::UpdateTrackParams,
 244         //       the current track parameters are changed
 245         //       and it could cause refit problems.
 246         //       We need to update only the following track parameters:
 247         //            the track length and expected times.
 248         //       Removed AliESDtrack::UpdateTrackParams call
 249         //       Called AliESDtrack::SetIntegratedTimes(...) and
 250         //       AliESDtrack::SetIntegratedLength() routines.
 251         /*
 252           AliTOFtrack *track = new AliTOFtrack(*seed);
 253           t->UpdateTrackParams(track,AliESDtrack::kTOFout); // to be checked - AdC
 254           delete track;
 255           Double_t time[10]; t->GetIntegratedTimes(time);
 256         */
 257
 258         Double_t time[10]; seed->GetIntegratedTimes(time);
 259         t->SetIntegratedTimes(time);
 260
 261         Double_t length =  seed->GetIntegratedLength();
 262         t->SetIntegratedLength(length);
 263
 264         Double_t alphaB = (Double_t)t->GetAlpha();
 265         Double_t xB = (Double_t)t->GetX();
 266         Double_t yB = (Double_t)t->GetY();
 267         Double_t zB = (Double_t)t->GetZ();
 268         Double_t p1B = (Double_t)t->GetSnp();
 269         Double_t p2B = (Double_t)t->GetTgl();
 270         Double_t p3B = (Double_t)t->GetSigned1Pt();
 271         const Double_t *covB = (Double_t*)t->GetCovariance();
 272         AliDebug(2,"Track params -now(before)-:");
 273         AliDebug(2,Form("    X: %f(%f), Y: %f(%f), Z: %f(%f) --- alpha: %f(%f)",
 274                         xB,xA,
 275                         yB,yA,
 276                         zB,zA,
 277                         alphaB,alphaA));
 278         AliDebug(2,Form("    p1: %f(%f), p2: %f(%f), p3: %f(%f)",
 279                         p1B,p1A,
 280                         p2B,p2A,
 281                         p3B,p3A));
 282         AliDebug(2,Form("    cov1: %f(%f), cov2: %f(%f), cov3: %f(%f)"
 283                         " cov4: %f(%f), cov5: %f(%f), cov6: %f(%f)"
 284                         " cov7: %f(%f), cov8: %f(%f), cov9: %f(%f)"
 285                         " cov10: %f(%f), cov11: %f(%f), cov12: %f(%f)"
 286                         " cov13: %f(%f), cov14: %f(%f), cov15: %f(%f)",
 287                         covB[0],covA[0],
 288                         covB[1],covA[1],
 289                         covB[2],covA[2],
 290                         covB[3],covA[3],
 291                         covB[4],covA[4],
 292                         covB[5],covA[5],
 293                         covB[6],covA[6],
 294                         covB[7],covA[7],
 295                         covB[8],covA[8],
 296                         covB[9],covA[9],
 297                         covB[10],covA[10],
 298                         covB[11],covA[11],
 299                         covB[12],covA[12],
 300                         covB[13],covA[13],
 301                         covB[14],covA[14]
 302                         ));
 303         AliDebug(3,Form(" TOF params: %6d  %f %f %f %f %f %6d %3d %f  %f %f %f %f %f",
 304                         i,
 305                         t->GetTOFsignalRaw(),
 306                         t->GetTOFsignal(),
 307                         t->GetTOFsignalToT(),
 308                         t->GetTOFsignalDz(),
 309                         t->GetTOFsignalDx(),
 310                         t->GetTOFCalChannel(),
 311                         t->GetTOFcluster(),
 312                         t->GetIntegratedLength(),
 313                         time[0], time[1], time[2], time[3], time[4]
 314                         )
 315                  );
 316       }
 317     }
 318   }
 319
 320   //Make TOF PID
 321   // Now done in AliESDpid
 322   // fPid->MakePID(event,timeZero);
 323
 324   fSeeds->Clear();
 325   //fTracks->Delete();
 326   fTracks->Clear();
 327   return 0;
 328
 329 }
 330 //_________________________________________________________________________
 331 void AliTOFtracker::CollectESD() {
 332    //prepare the set of ESD tracks to be matched to clusters in TOF
 333
 334   Int_t seedsTOF1=0;
 335   Int_t seedsTOF2=0;
 336
 337   TClonesArray &aTOFTrack = *fTracks;
 338   for (Int_t i=0; i<fNseeds; i++) {
 339
 340     AliESDtrack *t =(AliESDtrack*)fSeeds->At(i);
 341     if ((t->GetStatus()&AliESDtrack::kTPCout)==0)continue;
 342
 343     AliTOFtrack *track = new AliTOFtrack(*t); // New
 344     Float_t x = (Float_t)track->GetX(); //New
 345
 346     // TRD 'good' tracks, already propagated at 371 cm
 347     if ( ( (t->GetStatus()&AliESDtrack::kTRDout)!=0 ) &&
 348          ( x >= AliTOFGeometry::Rmin() ) ) {
 349       if  ( track->PropagateToInnerTOF() ) {
 350
 351         AliDebug(1,Form(" TRD propagated track till rho = %fcm."
 352                         " And then the track has been propagated till rho = %fcm.",
 353                         x, (Float_t)track->GetX()));
 354
 355         track->SetSeedIndex(i);
 356         t->UpdateTrackParams(track,AliESDtrack::kTOFin);
 357         new(aTOFTrack[fNseedsTOF]) AliTOFtrack(*track);
 358         fNseedsTOF++;
 359         seedsTOF1++;
 360       }
 361       delete track;
 362     }
 363
 364     // Propagate the rest of TPCbp
 365     else {
 366       if ( track->PropagateToInnerTOF() ) {
 367
 368         AliDebug(1,Form(" TPC propagated track till rho = %fcm."
 369                         " And then the track has been propagated till rho = %fcm.",
 370                         x, (Float_t)track->GetX()));
 371
 372         track->SetSeedIndex(i);
 373         t->UpdateTrackParams(track,AliESDtrack::kTOFin);
 374         new(aTOFTrack[fNseedsTOF]) AliTOFtrack(*track);
 375         fNseedsTOF++;
 376         seedsTOF2++;
 377       }
 378       delete track;
 379     }
 380   }
 381
 382   AliInfo(Form("Number of TOF seeds = %d (Type 1 = %d, Type 2 = %d)",fNseedsTOF,seedsTOF1,seedsTOF2));
 383
 384   // Sort according uncertainties on track position
 385   fTracks->Sort();
 386
 387 }
 388 //_________________________________________________________________________
 389 void AliTOFtracker::MatchTracks( Bool_t mLastStep){
 390
 391
 392   // Parameters used/regulating the reconstruction
 393
 394   static Float_t corrLen=0.75;
 395   static Float_t detDepth=15.3;
 396   static Float_t padDepth=0.5;
 397
 398   const Float_t kSpeedOfLight= 2.99792458e-2; // speed of light [cm/ps]
 399   const Float_t kTimeOffset = 0.; // time offset for tracking algorithm [ps]
 400
 401   Float_t dY=AliTOFGeometry::XPad();
 402   Float_t dZ=AliTOFGeometry::ZPad();
 403
 404   Float_t sensRadius = fkRecoParam->GetSensRadius();
 405   Float_t stepSize   = fkRecoParam->GetStepSize();
 406   Float_t scaleFact   = fkRecoParam->GetWindowScaleFact();
 407   Float_t dyMax=fkRecoParam->GetWindowSizeMaxY();
 408   Float_t dzMax=fkRecoParam->GetWindowSizeMaxZ();
 409   Float_t dCut=fkRecoParam->GetDistanceCut();
 410   Double_t maxChi2=fkRecoParam->GetMaxChi2TRD();
 411   Bool_t timeWalkCorr    = fkRecoParam->GetTimeWalkCorr();
 412   if(!mLastStep){
 413     AliDebug(1,"++++++++++++++TOF Reconstruction Parameters:++++++++++++ \n");
 414     AliDebug(1,Form("TOF sens radius: %f",sensRadius));
 415     AliDebug(1,Form("TOF step size: %f",stepSize));
 416     AliDebug(1,Form("TOF Window scale factor: %f",scaleFact));
 417     AliDebug(1,Form("TOF Window max dy: %f",dyMax));
 418     AliDebug(1,Form("TOF Window max dz: %f",dzMax));
 419     AliDebug(1,Form("TOF distance Cut: %f",dCut));
 420     AliDebug(1,Form("TOF Max Chi2: %f",maxChi2));
 421     AliDebug(1,Form("Time Walk Correction? : %d",timeWalkCorr));
 422   }
 423
 424   //Match ESD tracks to clusters in TOF
 425
 426   // Get the number of propagation steps
 427   Int_t nSteps=(Int_t)(detDepth/stepSize);
 428
 429   //PH Arrays (moved outside of the loop)
 430   Float_t * trackPos[4];
 431   for (Int_t ii=0; ii<4; ii++) trackPos[ii] = new Float_t[nSteps];
 432   Int_t * clind = new Int_t[fN];
 433
 434   // Some init
 435   const Int_t kNclusterMax = 1000; // related to fN value
 436   TGeoHMatrix global[kNclusterMax];
 437
 438   //The matching loop
 439   for (Int_t iseed=0; iseed<fNseedsTOF; iseed++) {
 440
 441     fTOFtrackPoints->Clear();
 442
 443     for (Int_t ii=0; ii<kNclusterMax; ii++)
 444       global[ii] = 0x0;
 445     AliTOFtrack *track =(AliTOFtrack*)fTracks->UncheckedAt(iseed);
 446     AliESDtrack *t =(AliESDtrack*)fSeeds->At(track->GetSeedIndex());
 447     //if ( t->GetTOFsignal()>0. ) continue;
 448     if ( (t->GetStatus()&AliESDtrack::kTOFout)!=0 ) continue;
 449     AliTOFtrack *trackTOFin = new AliTOFtrack(*track);
 450
 451     // Determine a window around the track
 452     Double_t x,par[5];
 453     trackTOFin->GetExternalParameters(x,par);
 454     Double_t cov[15];
 455     trackTOFin->GetExternalCovariance(cov);
 456
 457     if (cov[0]<0. || cov[2]<0.) {
 458       AliWarning(Form("Very strange track (%d)! At least one of its covariance matrix diagonal elements is negative!",iseed));
 459       //delete trackTOFin;
 460       //continue;
 461     }
 462
 463     Double_t dphi=
 464       scaleFact*
 465       ((5*TMath::Sqrt(TMath::Abs(cov[0])) + 0.5*dY + 2.5*TMath::Abs(par[2]))/sensRadius);
 466     Double_t dz=
 467        scaleFact*
 468        (5*TMath::Sqrt(TMath::Abs(cov[2])) + 0.5*dZ + 2.5*TMath::Abs(par[3]));
 469
 470     Double_t phi=TMath::ATan2(par[0],x) + trackTOFin->GetAlpha();
 471     if (phi<-TMath::Pi())phi+=2*TMath::Pi();
 472     if (phi>=TMath::Pi())phi-=2*TMath::Pi();
 473     Double_t z=par[1];
 474
 475     //upper limit on window's size.
 476     if (dz> dzMax) dz=dzMax;
 477     if (dphi*sensRadius> dyMax) dphi=dyMax/sensRadius;
 478
 479
 480     // find the clusters in the window of the track
 481     Int_t nc=0;
 482     for (Int_t k=FindClusterIndex(z-dz); k<fN; k++) {
 483
 484       if (nc>=kNclusterMax) {
 485         AliWarning("No more matchable clusters can be stored! Please, increase the corresponding vectors size.");
 486         break;
 487       }
 488
 489       AliTOFcluster *c=fClusters[k];
 490       if (c->GetZ() > z+dz) break;
 491       if (c->IsUsed()) continue;
 492       if (!c->GetStatus()) {
 493         AliDebug(1,"Cluster in channel declared bad!");
 494         continue; // skip bad channels as declared in OCDB
 495       }
 496
 497       Double_t dph=TMath::Abs(c->GetPhi()-phi);
 498       if (dph>TMath::Pi()) dph-=2.*TMath::Pi();
 499       if (TMath::Abs(dph)>dphi) continue;
 500
 501       Double_t yc=(c->GetPhi() - trackTOFin->GetAlpha())*c->GetR();
 502       Double_t p[2]={yc, c->GetZ()};
 503       Double_t cov2[3]= {dY*dY/12., 0., dZ*dZ/12.};
 504       if (trackTOFin->AliExternalTrackParam::GetPredictedChi2(p,cov2) > maxChi2)continue;
 505
 506       clind[nc] = k;
 507       Char_t path[200];
 508       Int_t ind[5];
 509       ind[0]=c->GetDetInd(0);
 510       ind[1]=c->GetDetInd(1);
 511       ind[2]=c->GetDetInd(2);
 512       ind[3]=c->GetDetInd(3);
 513       ind[4]=c->GetDetInd(4);
 514       fGeom->GetVolumePath(ind,path);
 515       gGeoManager->cd(path);
 516       global[nc] = *gGeoManager->GetCurrentMatrix();
 517       nc++;
 518     }
 519
 520     AliDebug(1,Form(" Number of matchable TOF clusters for the track number %d: %d",iseed,nc));
 521
 522     //start fine propagation
 523
 524     Int_t nStepsDone = 0;
 525     for( Int_t istep=0; istep<nSteps; istep++){
 526
 527       Float_t xs=AliTOFGeometry::RinTOF()+istep*0.1;
 528       Double_t ymax=xs*TMath::Tan(0.5*AliTOFGeometry::GetAlpha());
 529
 530       Bool_t skip=kFALSE;
 531       Double_t ysect=trackTOFin->GetYat(xs,skip);
 532       if (skip) break;
 533       if (ysect > ymax) {
 534         if (!trackTOFin->Rotate(AliTOFGeometry::GetAlpha())) {
 535           break;
 536         }
 537       } else if (ysect <-ymax) {
 538         if (!trackTOFin->Rotate(-AliTOFGeometry::GetAlpha())) {
 539           break;
 540         }
 541       }
 542
 543       if(!trackTOFin->PropagateTo(xs)) {
 544         break;
 545       }
 546
 547       nStepsDone++;
 548
 549       // store the running point (Globalrf) - fine propagation
 550
 551       Double_t r[3];
 552       trackTOFin->GetXYZ(r);
 553       trackPos[0][istep]= (Float_t) r[0];
 554       trackPos[1][istep]= (Float_t) r[1];
 555       trackPos[2][istep]= (Float_t) r[2];
 556       trackPos[3][istep]= trackTOFin->GetIntegratedLength();
 557     }
 558
 559
 560     Int_t nfound = 0;
 561     Bool_t accept = kFALSE;
 562     Bool_t isInside = kFALSE;
 563     for (Int_t istep=0; istep<nStepsDone; istep++) {
 564
 565       Float_t ctrackPos[3];
 566       ctrackPos[0] = trackPos[0][istep];
 567       ctrackPos[1] = trackPos[1][istep];
 568       ctrackPos[2] = trackPos[2][istep];
 569
 570       //now see whether the track matches any of the TOF clusters
 571
 572       Float_t dist3d[3];
 573       accept = kFALSE;
 574       for (Int_t i=0; i<nc; i++) {
 575         isInside = fGeom->IsInsideThePad((TGeoHMatrix*)(&global[i]),ctrackPos,dist3d);
 576
 577         if ( mLastStep ) {
 578           Float_t yLoc = dist3d[1];
 579           Float_t rLoc = TMath::Sqrt(dist3d[0]*dist3d[0]+dist3d[2]*dist3d[2]);
 580           accept = (TMath::Abs(yLoc)<padDepth*0.5 && rLoc<dCut);
 581           AliDebug(2," I am in the case mLastStep==kTRUE ");
 582         }
 583         else {
 584           accept = isInside;
 585         }
 586         if (accept) {
 587
 588           fTOFtrackPoints->AddLast(new AliTOFtrackPoint(clind[i],
 589                                                         TMath::Sqrt(dist3d[0]*dist3d[0] + dist3d[1]*dist3d[1] + dist3d[2]*dist3d[2]),
 590                                                         dist3d[2], dist3d[0],
 591                                                         AliTOFGeometry::RinTOF()+istep*0.1,trackPos[3][istep]));
 592
 593           AliDebug(2,Form(" dist3dLoc[0] = %f, dist3dLoc[1] = %f, dist3dLoc[2] = %f ",dist3d[0],dist3d[1],dist3d[2]));
 594           nfound++;
 595           if(accept &&!mLastStep)break;
 596         }//end if accept
 597
 598       } //end for on the clusters
 599       if(accept &&!mLastStep)break;
 600     } //end for on the steps
 601
 602     AliDebug(1,Form(" Number of track points for the track number %d: %d",iseed,nfound));
 603
 604
 605     if (nfound == 0 ) {
 606       fnunmatch++;
 607       delete trackTOFin;
 608       continue;
 609     }
 610
 611     fnmatch++;
 612
 613     // now choose the cluster to be matched with the track.
 614
 615     Int_t idclus=-1;
 616     Float_t  recL = 0.;
 617     Float_t  xpos=0.;
 618     Float_t  mindist=1000.;
 619     Float_t  mindistZ=0.;
 620     Float_t  mindistY=0.;
 621     for (Int_t iclus= 0; iclus<nfound;iclus++){
 622       AliTOFtrackPoint *matchableTOFcluster = (AliTOFtrackPoint*)fTOFtrackPoints->At(iclus);
 623       if (matchableTOFcluster->Distance()<mindist) {
 624         mindist = matchableTOFcluster->Distance();
 625         mindistZ = matchableTOFcluster->DistanceZ();  // Z distance in the
 626                                                       // RF of the hit pad
 627                                                       // closest to the
 628                                                       // reconstructed
 629                                                       // track
 630         mindistY = matchableTOFcluster->DistanceY(); // Y distance in the
 631                                                      // RF of the hit pad
 632                                                      // closest to the
 633                                                      // reconstructed
 634                                                      // track
 635         xpos = matchableTOFcluster->PropRadius();
 636         idclus = matchableTOFcluster->Index();
 637         recL = matchableTOFcluster->Length() + corrLen*0.5;
 638       }
 639     } // loop on found TOF track points
 640
 641     if (idclus==-1) {
 642       AliDebug(1,Form("Reconstructed track %d doesn't match any TOF cluster", iseed));
 643       delete trackTOFin;
 644       continue;
 645     }
 646
 647     AliTOFcluster *c=fClusters[idclus];
 648
 649     AliDebug(2, Form("%7d     %7d     %10d     %10d  %10d  %10d      %7d",
 650                      iseed,
 651                      fnmatch-1,
 652                      TMath::Abs(trackTOFin->GetLabel()),
 653                      c->GetLabel(0), c->GetLabel(1), c->GetLabel(2),
 654                      idclus)); // AdC
 655
 656     c->Use();
 657
 658     // Track length correction for matching Step 2
 659
 660     if (mLastStep) {
 661       Float_t rc = TMath::Sqrt(c->GetR()*c->GetR() + c->GetZ()*c->GetZ());
 662       Float_t rt = TMath::Sqrt(trackPos[0][70]*trackPos[0][70]
 663                                +trackPos[1][70]*trackPos[1][70]
 664                                +trackPos[2][70]*trackPos[2][70]);
 665       Float_t dlt=rc-rt;
 666       recL=trackPos[3][70]+dlt;
 667     }
 668
 669     if (
 670         (c->GetLabel(0)==TMath::Abs(trackTOFin->GetLabel()))
 671         ||
 672         (c->GetLabel(1)==TMath::Abs(trackTOFin->GetLabel()))
 673         ||
 674         (c->GetLabel(2)==TMath::Abs(trackTOFin->GetLabel()))
 675         ) {
 676       fngoodmatch++;
 677
 678        AliDebug(2,Form(" track label good %5d",trackTOFin->GetLabel()));
 679
 680     }
 681     else {
 682       fnbadmatch++;
 683
 684       AliDebug(2,Form(" track label  bad %5d",trackTOFin->GetLabel()));
 685
 686     }
 687
 688     delete trackTOFin;
 689
 690     //  Store quantities to be used in the TOF Calibration
 691     Float_t tToT=AliTOFGeometry::ToTBinWidth()*c->GetToT()*1E-3; // in ns
 692     t->SetTOFsignalToT(tToT);
 693     Float_t rawTime=AliTOFGeometry::TdcBinWidth()*c->GetTDCRAW()+kTimeOffset; // RAW time,in ps
 694     t->SetTOFsignalRaw(rawTime);
 695     t->SetTOFsignalDz(mindistZ);
 696     t->SetTOFsignalDx(mindistY);
 697     t->SetTOFDeltaBC(c->GetDeltaBC());
 698     t->SetTOFL0L1(c->GetL0L1Latency());
 699
 700     Float_t info[10] = {mindist,mindistY,mindistZ,
 701                         0.,0.,0.,0.,0.,0.,0.};
 702     t->SetTOFInfo(info);
 703     AliDebug(2,Form(" distance=%f; residual in the pad reference frame: dX=%f, dZ=%f", info[0],info[1],info[2]));
 704
 705
 706     Int_t ind[5];
 707     ind[0]=c->GetDetInd(0);
 708     ind[1]=c->GetDetInd(1);
 709     ind[2]=c->GetDetInd(2);
 710     ind[3]=c->GetDetInd(3);
 711     ind[4]=c->GetDetInd(4);
 712     Int_t calindex = AliTOFGeometry::GetIndex(ind);
 713     t->SetTOFCalChannel(calindex);
 714
 715     // keep track of the track labels in the matched cluster
 716     Int_t tlab[3];
 717     tlab[0]=c->GetLabel(0);
 718     tlab[1]=c->GetLabel(1);
 719     tlab[2]=c->GetLabel(2);
 720     AliDebug(2,Form(" tdc time of the matched track %6d = ",c->GetTDC()));
 721     Double_t tof=AliTOFGeometry::TdcBinWidth()*c->GetTDC()+kTimeOffset; // in ps
 722     AliDebug(2,Form(" tof time of the matched track: %f = ",tof));
 723     Double_t tofcorr=tof;
 724     if(timeWalkCorr)tofcorr=CorrectTimeWalk(mindistZ,tof);
 725     AliDebug(2,Form(" tof time of the matched track, after TW corr: %f = ",tofcorr));
 726     //Set TOF time signal and pointer to the matched cluster
 727     t->SetTOFsignal(tofcorr);
 728     t->SetTOFcluster(idclus); // pointing to the recPoints tree
 729
 730     AliDebug(2,Form(" Setting TOF raw time: %f, z distance: %f  corrected time: %f ",rawTime,mindistZ,tofcorr));
 731
 732     //Tracking info
 733     Double_t time[AliPID::kSPECIES]; t->GetIntegratedTimes(time); // in ps
 734     Double_t mom=t->GetP();
 735     AliDebug(2,Form(" Momentum for track %d -> %f", iseed,mom));
 736     for (Int_t j=0;j<AliPID::kSPECIES;j++) {
 737       Double_t mass=AliPID::ParticleMass(j);
 738       time[j]+=(recL-trackPos[3][0])/kSpeedOfLight*TMath::Sqrt(mom*mom+mass*mass)/mom;
 739     }
 740
 741     AliTOFtrack *trackTOFout = new AliTOFtrack(*t);
 742     trackTOFout->PropagateTo(xpos);
 743
 744     // Fill the track residual histograms.
 745     FillResiduals(trackTOFout,c,kFALSE);
 746
 747     t->UpdateTrackParams(trackTOFout,AliESDtrack::kTOFout);
 748     t->SetIntegratedLength(recL);
 749     t->SetIntegratedTimes(time);
 750     t->SetTOFLabel(tlab);
 751
 752
 753     // Fill Reco-QA histos for Reconstruction
 754     fHRecNClus->Fill(nc);
 755     fHRecDist->Fill(mindist);
 756     if (cov[0]>=0.)
 757       fHRecSigYVsP->Fill(mom,TMath::Sqrt(cov[0]));
 758     else
 759       fHRecSigYVsP->Fill(mom,-TMath::Sqrt(-cov[0]));
 760     if (cov[2]>=0.)
 761       fHRecSigZVsP->Fill(mom,TMath::Sqrt(cov[2]));
 762     else
 763       fHRecSigZVsP->Fill(mom,-TMath::Sqrt(-cov[2]));
 764     fHRecSigYVsPWin->Fill(mom,dphi*sensRadius);
 765     fHRecSigZVsPWin->Fill(mom,dz);
 766
 767     // Fill Tree for on-the-fly offline Calibration
 768
 769     if ( !((t->GetStatus() & AliESDtrack::kTIME)==0 ) ) {
 770       fIch=calindex;
 771       fToT=tToT;
 772       fTime=rawTime;
 773       fExpTimePi=time[2];
 774       fExpTimeKa=time[3];
 775       fExpTimePr=time[4];
 776       fCalTree->Fill();
 777     }
 778     delete trackTOFout;
 779   }
 780
 781   for (Int_t ii=0; ii<4; ii++) delete [] trackPos[ii];
 782   delete [] clind;
 783
 784 }
 785 //_________________________________________________________________________
 786 Int_t AliTOFtracker::LoadClusters(TTree *cTree) {
 787   //--------------------------------------------------------------------
 788   //This function loads the TOF clusters
 789   //--------------------------------------------------------------------
 790
 791   Int_t npadX = AliTOFGeometry::NpadX();
 792   Int_t npadZ = AliTOFGeometry::NpadZ();
 793   Int_t nStripA = AliTOFGeometry::NStripA();
 794   Int_t nStripB = AliTOFGeometry::NStripB();
 795   Int_t nStripC = AliTOFGeometry::NStripC();
 796
 797   TBranch *branch=cTree->GetBranch("TOF");
 798   if (!branch) {
 799     AliError("can't get the branch with the TOF clusters !");
 800     return 1;
 801   }
 802
 803   static TClonesArray dummy("AliTOFcluster",10000);
 804   dummy.Clear();
 805   TClonesArray *clusters=&dummy;
 806   branch->SetAddress(&clusters);
 807
 808   cTree->GetEvent(0);
 809   Int_t nc=clusters->GetEntriesFast();
 810   fHDigNClus->Fill(nc);
 811
 812   AliInfo(Form("Number of clusters: %d",nc));
 813
 814   for (Int_t i=0; i<nc; i++) {
 815     AliTOFcluster *c=(AliTOFcluster*)clusters->UncheckedAt(i);
 816 //PH    fClusters[i]=new AliTOFcluster(*c); fN++;
 817     fClusters[i]=c; fN++;
 818
 819   // Fill Digits QA histos
 820
 821     Int_t isector = c->GetDetInd(0);
 822     Int_t iplate = c->GetDetInd(1);
 823     Int_t istrip = c->GetDetInd(2);
 824     Int_t ipadX = c->GetDetInd(4);
 825     Int_t ipadZ = c->GetDetInd(3);
 826
 827     Float_t time =(AliTOFGeometry::TdcBinWidth()*c->GetTDC())*1E-3; // in ns
 828     Float_t tot = (AliTOFGeometry::TdcBinWidth()*c->GetToT())*1E-3;//in ns
 829
 830     Int_t stripOffset = 0;
 831     switch (iplate) {
 832     case 0:
 833       stripOffset = 0;
 834       break;
 835     case 1:
 836       stripOffset = nStripC;
 837       break;
 838     case 2:
 839       stripOffset = nStripC+nStripB;
 840       break;
 841     case 3:
 842       stripOffset = nStripC+nStripB+nStripA;
 843       break;
 844     case 4:
 845       stripOffset = nStripC+nStripB+nStripA+nStripB;
 846       break;
 847     default:
 848       AliError(Form("Wrong plate number in TOF (%d) !",iplate));
 849       break;
 850     };
 851     Int_t zindex=npadZ*(istrip+stripOffset)+(ipadZ+1);
 852     Int_t phiindex=npadX*isector+ipadX+1;
 853     fHDigClusMap->Fill(zindex,phiindex);
 854     fHDigClusTime->Fill(time);
 855     fHDigClusToT->Fill(tot);
 856
 857   }
 858
 859
 860   return 0;
 861 }
 862 //_________________________________________________________________________
 863 void AliTOFtracker::UnloadClusters() {
 864   //--------------------------------------------------------------------
 865   //This function unloads TOF clusters
 866   //--------------------------------------------------------------------
 867   for (Int_t i=0; i<fN; i++) {
 868 //PH    delete fClusters[i];
 869     fClusters[i] = 0x0;
 870   }
 871   fN=0;
 872 }
 873
 874 //_________________________________________________________________________
 875 Int_t AliTOFtracker::FindClusterIndex(Double_t z) const {
 876   //--------------------------------------------------------------------
 877   // This function returns the index of the nearest cluster
 878   //--------------------------------------------------------------------
 879   if (fN==0) return 0;
 880   if (z <= fClusters[0]->GetZ()) return 0;
 881   if (z > fClusters[fN-1]->GetZ()) return fN;
 882   Int_t b=0, e=fN-1, m=(b+e)/2;
 883   for (; b<e; m=(b+e)/2) {
 884     if (z > fClusters[m]->GetZ()) b=m+1;
 885     else e=m;
 886   }
 887   return m;
 888 }
 889
 890 //_________________________________________________________________________
 891 Bool_t AliTOFtracker::GetTrackPoint(Int_t index, AliTrackPoint& p) const
 892 {
 893   // Get track space point with index i
 894   // Coordinates are in the global system
 895   AliTOFcluster *cl = fClusters[index];
 896   Float_t xyz[3];
 897   xyz[0] = cl->GetR()*TMath::Cos(cl->GetPhi());
 898   xyz[1] = cl->GetR()*TMath::Sin(cl->GetPhi());
 899   xyz[2] = cl->GetZ();
 900   Float_t phiangle = (Int_t(cl->GetPhi()*TMath::RadToDeg()/20.)+0.5)*20.*TMath::DegToRad();
 901   Float_t sinphi = TMath::Sin(phiangle), cosphi = TMath::Cos(phiangle);
 902   Float_t tiltangle = AliTOFGeometry::GetAngles(cl->GetDetInd(1),cl->GetDetInd(2))*TMath::DegToRad();
 903   Float_t sinth = TMath::Sin(tiltangle), costh = TMath::Cos(tiltangle);
 904   Float_t sigmay2 = AliTOFGeometry::XPad()*AliTOFGeometry::XPad()/12.;
 905   Float_t sigmaz2 = AliTOFGeometry::ZPad()*AliTOFGeometry::ZPad()/12.;
 906   Float_t cov[6];
 907   cov[0] = sinphi*sinphi*sigmay2 + cosphi*cosphi*sinth*sinth*sigmaz2;
 908   cov[1] = -sinphi*cosphi*sigmay2 + sinphi*cosphi*sinth*sinth*sigmaz2;
 909   cov[2] = -cosphi*sinth*costh*sigmaz2;
 910   cov[3] = cosphi*cosphi*sigmay2 + sinphi*sinphi*sinth*sinth*sigmaz2;
 911   cov[4] = -sinphi*sinth*costh*sigmaz2;
 912   cov[5] = costh*costh*sigmaz2;
 913   p.SetXYZ(xyz[0],xyz[1],xyz[2],cov);
 914
 915   // Detector numbering scheme
 916   Int_t nSector = AliTOFGeometry::NSectors();
 917   Int_t nPlate  = AliTOFGeometry::NPlates();
 918   Int_t nStripA = AliTOFGeometry::NStripA();
 919   Int_t nStripB = AliTOFGeometry::NStripB();
 920   Int_t nStripC = AliTOFGeometry::NStripC();
 921
 922   Int_t isector = cl->GetDetInd(0);
 923   if (isector >= nSector)
 924     AliError(Form("Wrong sector number in TOF (%d) !",isector));
 925   Int_t iplate = cl->GetDetInd(1);
 926   if (iplate >= nPlate)
 927     AliError(Form("Wrong plate number in TOF (%d) !",iplate));
 928   Int_t istrip = cl->GetDetInd(2);
 929
 930   Int_t stripOffset = 0;
 931   switch (iplate) {
 932   case 0:
 933     stripOffset = 0;
 934     break;
 935   case 1:
 936     stripOffset = nStripC;
 937     break;
 938   case 2:
 939     stripOffset = nStripC+nStripB;
 940     break;
 941   case 3:
 942     stripOffset = nStripC+nStripB+nStripA;
 943     break;
 944   case 4:
 945     stripOffset = nStripC+nStripB+nStripA+nStripB;
 946     break;
 947   default:
 948     AliError(Form("Wrong plate number in TOF (%d) !",iplate));
 949     break;
 950   };
 951
 952   Int_t idet = (2*(nStripC+nStripB)+nStripA)*isector +
 953                stripOffset +
 954                istrip;
 955   UShort_t volid = AliGeomManager::LayerToVolUID(AliGeomManager::kTOF,idet);
 956   p.SetVolumeID((UShort_t)volid);
 957   return kTRUE;
 958 }
 959 //_________________________________________________________________________
 960 void AliTOFtracker::InitCheckHists() {
 961
 962   //Init histos for Digits/Reco QA and Calibration
 963
 964
 965   TDirectory *dir = gDirectory;
 966   TFile *logFileTOF = 0;
 967
 968   TSeqCollection *list = gROOT->GetListOfFiles();
 969   int n = list->GetEntries();
 970   Bool_t isThere=kFALSE;
 971   for(int i=0; i<n; i++) {
 972     logFileTOF = (TFile*)list->At(i);
 973     if (strstr(logFileTOF->GetName(), "TOFQA.root")){
 974       isThere=kTRUE;
 975       break;
 976     }
 977   }
 978
 979   if(!isThere)logFileTOF = new TFile( "TOFQA.root","RECREATE");
 980   logFileTOF->cd();
 981
 982   fCalTree = new TTree("CalTree", "Tree for TOF calibration");
 983   fCalTree->Branch("TOFchannelindex",&fIch,"iTOFch/I");
 984   fCalTree->Branch("ToT",&fToT,"TOFToT/F");
 985   fCalTree->Branch("TOFtime",&fTime,"TOFtime/F");
 986   fCalTree->Branch("PionExpTime",&fExpTimePi,"PiExpTime/F");
 987   fCalTree->Branch("KaonExpTime",&fExpTimeKa,"KaExpTime/F");
 988   fCalTree->Branch("ProtonExpTime",&fExpTimePr,"PrExpTime/F");
 989
 990   //Digits "QA"
 991   fHDigClusMap = new TH2F("TOFDig_ClusMap", "",182,0.5,182.5,864, 0.5,864.5);
 992   fHDigNClus = new TH1F("TOFDig_NClus", "",200,0.5,200.5);
 993   fHDigClusTime = new TH1F("TOFDig_ClusTime", "",2000,0.,200.);
 994   fHDigClusToT = new TH1F("TOFDig_ClusToT", "",500,0.,100);
 995
 996   //Reco "QA"
 997   fHRecNClus =new TH1F("TOFRec_NClusW", "",50,0.5,50.5);
 998   fHRecDist=new TH1F("TOFRec_Dist", "",50,0.5,10.5);
 999   fHRecSigYVsP=new TH2F("TOFDig_SigYVsP", "",40,0.,4.,100, 0.,5.);
1000   fHRecSigZVsP=new TH2F("TOFDig_SigZVsP", "",40,0.,4.,100, 0.,5.);
1001   fHRecSigYVsPWin=new TH2F("TOFDig_SigYVsPWin", "",40,0.,4.,100, 0.,50.);
1002   fHRecSigZVsPWin=new TH2F("TOFDig_SigZVsPWin", "",40,0.,4.,100, 0.,50.);
1003
1004   dir->cd();
1005
1006 }
1007
1008 //_________________________________________________________________________
1009 void AliTOFtracker::SaveCheckHists() {
1010
1011   //write histos for Digits/Reco QA and Calibration
1012
1013   TDirectory *dir = gDirectory;
1014   TFile *logFileTOF = 0;
1015
1016   TSeqCollection *list = gROOT->GetListOfFiles();
1017   int n = list->GetEntries();
1018   Bool_t isThere=kFALSE;
1019   for(int i=0; i<n; i++) {
1020     logFileTOF = (TFile*)list->At(i);
1021     if (strstr(logFileTOF->GetName(), "TOFQA.root")){
1022       isThere=kTRUE;
1023       break;
1024     }
1025   }
1026
1027   if(!isThere) {
1028           AliError(Form("File TOFQA.root not found!! not wring histograms...."));
1029           return;
1030   }
1031   logFileTOF->cd();
1032   fHDigClusMap->Write(fHDigClusMap->GetName(), TObject::kOverwrite);
1033   fHDigNClus->Write(fHDigNClus->GetName(), TObject::kOverwrite);
1034   fHDigClusTime->Write(fHDigClusTime->GetName(), TObject::kOverwrite);
1035   fHDigClusToT->Write(fHDigClusToT->GetName(), TObject::kOverwrite);
1036   fHRecNClus->Write(fHRecNClus->GetName(), TObject::kOverwrite);
1037   fHRecDist->Write(fHRecDist->GetName(), TObject::kOverwrite);
1038   fHRecSigYVsP->Write(fHRecSigYVsP->GetName(), TObject::kOverwrite);
1039   fHRecSigZVsP->Write(fHRecSigZVsP->GetName(), TObject::kOverwrite);
1040   fHRecSigYVsPWin->Write(fHRecSigYVsPWin->GetName(), TObject::kOverwrite);
1041   fHRecSigZVsPWin->Write(fHRecSigZVsPWin->GetName(), TObject::kOverwrite);
1042   fCalTree->Write(fCalTree->GetName(),TObject::kOverwrite);
1043   logFileTOF->Flush();
1044
1045   dir->cd();
1046   }
1047 //_________________________________________________________________________
1048 Float_t AliTOFtracker::CorrectTimeWalk( Float_t dist, Float_t tof) const {
1049
1050   //dummy, for the moment
1051   Float_t tofcorr=0.;
1052   if(dist<AliTOFGeometry::ZPad()*0.5){
1053     tofcorr=tof;
1054     //place here the actual correction
1055   }else{
1056     tofcorr=tof;
1057   }
1058   return tofcorr;
1059 }
1060 //_________________________________________________________________________
1061
1062 void AliTOFtracker::FillClusterArray(TObjArray* arr) const
1063 {
1064   //
1065   // Returns the TOF cluster array
1066   //
1067
1068   if (fN==0)
1069     arr = 0x0;
1070   else
1071     for (Int_t i=0; i<fN; ++i) arr->Add(fClusters[i]);
1072
1073 }
1074 //_________________________________________________________________________
1075