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