MFT/AliMuonForwardTrackFinder.cxx

   1 // ROOT includes
   2 #include "TObject.h"
   3 #include "TClonesArray.h"
   4 #include "TObjArray.h"
   5 #include "TH2D.h"
   6 #include "TH1D.h"
   7 #include "TFile.h"
   8 #include "TGeoManager.h"
   9 #include "TMatrixD.h"
  10 #include "TParticle.h"
  11 #include "TMath.h"
  12 #include "TGraph.h"
  13 #include "TEllipse.h"
  14 #include "TCanvas.h"
  15 #include "TString.h"
  16 #include "TLatex.h"
  17 #include "TMarker.h"
  18 #include "TNtuple.h"
  19 #include "TRandom.h"
  20 #include "TIterator.h"
  21
  22 // STEER includes
  23 #include "AliLog.h"
  24 #include "AliRun.h"
  25 #include "AliRunLoader.h"
  26 #include "AliLoader.h"
  27 #include "AliHeader.h"
  28 #include "AliMC.h"
  29 #include "AliStack.h"
  30 #include "AliMagF.h"
  31 #include "AliTracker.h"
  32 #include "AliGRPObject.h"
  33 #include "AliCDBEntry.h"
  34 #include "AliCDBManager.h"
  35
  36 // MUON includes
  37 #include "AliMUONConstants.h"
  38 #include "AliMUONTrack.h"
  39 #include "AliMUONRecoCheck.h"
  40 #include "AliMUONTrackParam.h"
  41 #include "AliMUONTrackExtrap.h"
  42 #include "AliMUONVTrackStore.h"
  43 #include "AliMUONVCluster.h"
  44
  45 // MFT includes
  46 #include "AliMuonForwardTrack.h"
  47 #include "AliMFTCluster.h"
  48 #include "AliMFT.h"
  49 #include "AliMFTSegmentation.h"
  50 #include "AliMFTConstants.h"
  51
  52 #include "AliMuonForwardTrackFinder.h"
  53
  54 //====================================================================================================================================================
  55 //
  56 // Class for the creation of the "global muon tracks" built from the clusters in the
  57 // muon spectrometer and the clusters of the Muon Forward Tracker. QA histograms are also created
  58 //
  59 // Contact author: antonio.uras@cern.ch
  60 //
  61 //====================================================================================================================================================
  62
  63 const Double_t AliMuonForwardTrackFinder::fRadLengthSi = AliMFTConstants::fRadLengthSi;
  64
  65 ClassImp(AliMuonForwardTrackFinder)
  66
  67 //=====================================================================================================
  68
  69 AliMuonForwardTrackFinder::AliMuonForwardTrackFinder():
  70   TObject(),
  71   fRun(0),
  72   fNEventsToAnalyze(0),
  73   fSigmaClusterCut(0),
  74   fChi2GlobalCut(0),
  75   fSigmaSpectrometerCut(0),
  76   fVertexErrorX(0.015),
  77   fVertexErrorY(0.015),
  78   fVertexErrorZ(0.010),
  79   fNFinalCandidatesCut(0),
  80   fReadDir(0),
  81   fOutDir(0),
  82   fDrawOption(0),
  83
  84   fDistanceFromGoodClusterAndTrackAtLastPlane(-1),
  85   fDistanceFromBestClusterAndTrackAtLastPlane(-1),
  86
  87   fRAbsorberCut(0),
  88   fLowPtCut(0),
  89   fNPlanesMFT(0),
  90   fNPlanesMFTAnalyzed(0),
  91   fNMaxMissingMFTClusters(0),
  92
  93   fEv(0),
  94   fLabelMC(0),
  95
  96   fHistRadiusEndOfAbsorber(0),
  97   fHistNGoodClustersForFinalTracks(0),
  98   fHistDistanceGoodClusterFromTrackMinusDistanceBestClusterFromTrackAtLastPlane(0),
  99   fHistDistanceGoodClusterFromTrackAtLastPlane(0),
 100
 101   fNtuFinalCandidates(0),
 102   fNtuFinalBestCandidates(0),
 103
 104   fCanvas(0),
 105
 106   fTxtMuonHistory(0),
 107   fTxtTrackGoodClusters(0),
 108   fTxtTrackFinalChi2(0),
 109   fTxtTrackMomentum(0),
 110   fTxtFinalCandidates(0),
 111   fTxtDummy(0),
 112   fTxtAllClust(0),
 113   fTxtClustGoodChi2(0),
 114   fTxtClustMC(0),
 115   fTxtClustOfTrack(0),
 116   fMrkAllClust(0),
 117   fMrkClustGoodChi2(0),
 118   fMrkClustMC(0),
 119   fMrkClustOfTrack(0),
 120
 121   fCountRealTracksAnalyzed(0),
 122   fMaxNTracksToBeAnalyzed(99999999),
 123   fCountRealTracksWithRefMC(0),
 124   fCountRealTracksWithRefMC_andTrigger(0),
 125   fCountRealTracksWithRefMC_andTrigger_andGoodPt(0),
 126   fCountRealTracksWithRefMC_andTrigger_andGoodPt_andGoodTheta(0),
 127   fCountRealTracksAnalyzedOfEvent(0),
 128   fCountRealTracksAnalyzedWithFinalCandidates(0),
 129
 130   fFileCluster(0),
 131   fFileESD(0),
 132   fFile_gAlice(0),
 133
 134   fRunLoader(0),
 135   fMFTLoader(0),
 136   fMuonRecoCheck(0),
 137   fMFTClusterTree(0),
 138   fMuonTrackReco(0),
 139   fCurrentTrack(0),
 140   fFinalBestCandidate(0),
 141   fIsCurrentMuonTrackable(0),
 142   fCandidateTracks(0),
 143   fTrackStore(0),
 144   fTrackRefStore(0),
 145   fNextTrack(0),
 146   fStack(0),
 147   fMFT(0),
 148   fSegmentation(0),
 149   fOutputTreeFile(0),
 150   fOutputQAFile(0),
 151   fOutputEventTree(0),
 152   fMuonForwardTracks(0),
 153   fMatchingMode(-1),
 154   fGRPData(0),
 155   fRunInfo(0),
 156   fBransonCorrection(kTRUE)
 157
 158 {
 159
 160   // Default constructor
 161
 162   for (Int_t iPlane=0; iPlane<AliMFTConstants::fNMaxPlanes; iPlane++) {
 163
 164     fHistNTracksAfterExtrapolation[iPlane] = 0;
 165     fHistChi2Cluster_GoodCluster[iPlane] = 0;
 166     fHistChi2Cluster_BadCluster[iPlane] = 0;
 167     fHistResearchRadius[iPlane] = 0;
 168
 169     fIsGoodClusterInPlane[iPlane] = kFALSE;
 170
 171     fHistChi2Cluster_GoodCluster[iPlane] = 0;
 172     fHistChi2Cluster_BadCluster[iPlane]  = 0;
 173
 174     fHistGlobalChi2AtPlaneFor_GOOD_CandidatesOfTrackableMuons[iPlane] = 0;
 175     fHistGlobalChi2AtPlaneFor_BAD_CandidatesOfTrackableMuons[iPlane] = 0;
 176
 177     fZPlane[iPlane] = 0.;
 178     fRPlaneMax[iPlane] = 0.;
 179     fRPlaneMin[iPlane] = 0.;
 180
 181     for (Int_t i=0; i<4; i++) fGrMFTPlane[i][iPlane] = 0;
 182     fCircleExt[iPlane] = 0;
 183     fCircleInt[iPlane] = 0;
 184
 185     fTxtTrackChi2[iPlane] = 0;
 186
 187     fIsClusterCompatible[iPlane] = 0;
 188
 189     fMFTClusterArray[iPlane]      = 0;
 190     fMFTClusterArrayFront[iPlane] = new TClonesArray("AliMFTCluster");
 191     fMFTClusterArrayBack[iPlane]  = new TClonesArray("AliMFTCluster");
 192
 193     fIsPlaneMandatory[iPlane] = kFALSE;
 194
 195     fMinResearchRadiusAtPlane[iPlane] = 0.;
 196
 197   }
 198
 199   //  fNextTrack = 0;
 200
 201   fHistDistanceGoodClusterFromTrackMinusDistanceBestClusterFromTrackAtLastPlane = 0;
 202   fHistDistanceGoodClusterFromTrackAtLastPlane = 0;
 203
 204   fMFTClusterTree = 0;
 205   fCandidateTracks = 0;
 206
 207   fOutputTreeFile    = new TFile("MuonGlobalTracks.root", "recreate");
 208   fOutputEventTree   = new TTree("AliMuonForwardTracks", "Tree of AliMuonForwardTracks");
 209   fMuonForwardTracks = new TClonesArray("AliMuonForwardTrack");
 210   fOutputEventTree   -> Branch("tracks", &fMuonForwardTracks);
 211
 212 }
 213
 214 //=====================================================================================================
 215
 216 AliMuonForwardTrackFinder::~AliMuonForwardTrackFinder() {
 217
 218   for (Int_t iPlane=0; iPlane<AliMFTConstants::fNMaxPlanes; iPlane++) {
 219
 220     delete fHistNTracksAfterExtrapolation[iPlane];
 221     delete fHistChi2Cluster_GoodCluster[iPlane];
 222     delete fHistChi2Cluster_BadCluster[iPlane];
 223     delete fHistResearchRadius[iPlane];
 224
 225     delete fHistChi2Cluster_GoodCluster[iPlane];
 226     delete fHistChi2Cluster_BadCluster[iPlane];
 227
 228     delete fHistGlobalChi2AtPlaneFor_GOOD_CandidatesOfTrackableMuons[iPlane];
 229     delete fHistGlobalChi2AtPlaneFor_BAD_CandidatesOfTrackableMuons[iPlane];
 230
 231     for (Int_t i=0; i<4; i++) delete fGrMFTPlane[i][iPlane];
 232     delete fCircleExt[iPlane];
 233     delete fCircleInt[iPlane];
 234
 235     delete fTxtTrackChi2[iPlane];
 236
 237     delete fMFTClusterArray[iPlane];
 238     delete fMFTClusterArrayFront[iPlane];
 239     delete fMFTClusterArrayBack[iPlane];
 240
 241   }
 242
 243   delete fNtuFinalCandidates;
 244   delete fNtuFinalBestCandidates;
 245
 246   delete fHistRadiusEndOfAbsorber;
 247
 248   delete fHistNGoodClustersForFinalTracks;
 249   delete fHistDistanceGoodClusterFromTrackMinusDistanceBestClusterFromTrackAtLastPlane;         //
 250   delete fHistDistanceGoodClusterFromTrackAtLastPlane;                                          //
 251
 252   delete fCanvas;
 253
 254   delete fTxtMuonHistory;
 255   delete fTxtTrackGoodClusters;
 256   delete fTxtTrackFinalChi2;
 257   delete fTxtTrackMomentum;
 258   delete fTxtFinalCandidates;
 259   delete fTxtDummy;
 260   delete fTxtAllClust;
 261   delete fTxtClustGoodChi2;
 262   delete fTxtClustMC;
 263   delete fTxtClustOfTrack;
 264   delete fMrkAllClust;
 265   delete fMrkClustGoodChi2;
 266   delete fMrkClustMC;
 267   delete fMrkClustOfTrack;
 268
 269   delete fFileCluster;
 270   delete fFileESD;
 271   delete fFile_gAlice;
 272
 273   delete fRunLoader;
 274   delete fMFTLoader;
 275   delete fMuonRecoCheck;
 276
 277   delete fMFTClusterTree;
 278
 279   delete fMuonTrackReco;
 280   delete fCurrentTrack;
 281   delete fFinalBestCandidate;
 282
 283   delete fCandidateTracks;
 284
 285   delete fTrackStore;
 286   delete fTrackRefStore;
 287
 288   delete fNextTrack;
 289
 290   delete fStack;
 291
 292   delete fMFT;
 293   delete fSegmentation;
 294
 295   delete fOutputTreeFile;
 296   delete fOutputQAFile;
 297   delete fOutputEventTree;
 298
 299   delete fMuonForwardTracks;
 300
 301   delete fGRPData;
 302   delete fRunInfo;
 303
 304 }
 305
 306 //=====================================================================================================
 307
 308 void AliMuonForwardTrackFinder::Init(Int_t nRun,
 309                                      Char_t *readDir,
 310                                      Char_t *outDir,
 311                                      Int_t nEventsToAnalyze) {
 312
 313   if (fRunLoader) {
 314     AliInfo("WARNING: run already initialized!!\n");
 315   }
 316
 317   SetRun(nRun);
 318   SetReadDir(readDir);
 319   SetOutDir(outDir);
 320
 321   AliInfo(Form("input  dir = %s\n", fReadDir.Data()));
 322   AliInfo(Form("output dir = %s\n", fOutDir.Data()));
 323
 324   // -------------------------- initializing files...
 325
 326   AliInfo(Form("initializing files for run %d...\n", fRun));
 327
 328   Char_t geoFileName[300];
 329   Char_t esdFileName[300];
 330   Char_t gAliceName[300];
 331   Char_t clusterName[300];
 332
 333   snprintf(geoFileName , 300, "%s/geometry.root",      fReadDir.Data());
 334   snprintf(esdFileName , 300, "%s/AliESDs.root" ,      fReadDir.Data());
 335   snprintf(gAliceName  , 300, "%s/galice.root"  ,      fReadDir.Data());
 336   snprintf(clusterName , 300, "%s/MFT.RecPoints.root", fReadDir.Data());
 337
 338   // Import TGeo geometry (needed by AliMUONTrackExtrap::ExtrapToVertex)
 339   if (!gGeoManager) {
 340     TGeoManager::Import(geoFileName);
 341     if (!gGeoManager) {
 342       AliError(Form("getting geometry from file %s failed", geoFileName));
 343       return;
 344     }
 345   }
 346
 347   fFileESD = new TFile(esdFileName);
 348   if (!fFileESD || !fFileESD->IsOpen()) return;
 349   else AliInfo(Form("file %s successfully opened\n", fFileESD->GetName()));
 350
 351   fMuonRecoCheck = new AliMUONRecoCheck(esdFileName, Form("%s/generated/", fReadDir.Data()));       // Utility class to check reconstruction
 352   fFile_gAlice = new TFile(gAliceName);
 353   if (!fFile_gAlice || !fFile_gAlice->IsOpen()) return;
 354   else AliInfo(Form("file %s successfully opened\n", fFile_gAlice->GetName()));
 355
 356   fRunLoader = AliRunLoader::Open(gAliceName);
 357   gAlice = fRunLoader->GetAliRun();
 358   if (!gAlice) fRunLoader->LoadgAlice();
 359   fMFT = (AliMFT*) gAlice->GetDetector("MFT");
 360   fSegmentation = fMFT->GetSegmentation();
 361   SetNPlanesMFT(fSegmentation->GetNPlanes());
 362
 363   if (!SetRunNumber()) return;
 364   if (!InitGRP()) return;
 365   AliMUONTrackExtrap::SetField();        // set the magnetic field for track extrapolations
 366
 367   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
 368     fZPlane[iPlane]    = fSegmentation->GetPlane(iPlane)->GetZCenter();
 369     fRPlaneMax[iPlane] = fSegmentation->GetPlane(iPlane)->GetRMaxSupport();
 370     fRPlaneMin[iPlane] = fSegmentation->GetPlane(iPlane)->GetRMinSupport();
 371   }
 372
 373   // Loading MFT clusters
 374   fMFTLoader = fRunLoader->GetDetectorLoader("MFT");
 375   fMFTLoader->LoadRecPoints("READ");
 376
 377   fMFTClusterTree = fMFTLoader->TreeR();
 378
 379   Int_t nEventsInFile = fMuonRecoCheck->NumberOfEvents();
 380   if (!nEventsInFile) {
 381     AliError("no events available!!!\n");
 382     return;
 383   }
 384   if (nEventsInFile<nEventsToAnalyze || nEventsToAnalyze<0) fNEventsToAnalyze = nEventsInFile;
 385   else fNEventsToAnalyze = nEventsToAnalyze;
 386
 387   fCandidateTracks = new TClonesArray("AliMuonForwardTrack",50000);
 388
 389   // -------------------------- initializing histograms...
 390
 391   AliInfo("\ninitializing histograms...\n");
 392   BookHistos();
 393   SetTitleHistos();
 394   AliInfo("... done!\n\n");
 395
 396   // -------------------------- initializing graphics...
 397
 398   AliInfo("initializing graphics...\n");
 399   BookPlanes();
 400   AliInfo("... done!\n\n");
 401
 402   SetSigmaSpectrometerCut(4.0);
 403   SetSigmaClusterCut(4.5);
 404   SetChi2GlobalCut(2.0);
 405   SetNFinalCandidatesCut(10);
 406   SetRAbsorberCut(26.4);
 407   SetLowPtCut(0.5);
 408
 409 }
 410
 411 //======================================================================================================================================
 412
 413 Bool_t AliMuonForwardTrackFinder::LoadNextEvent() {
 414
 415   // load next reconstructed event from the tree
 416
 417   if (fEv) FillOutputTree();
 418
 419   if (fEv>=fNEventsToAnalyze) return kFALSE;
 420
 421   fCountRealTracksAnalyzedOfEvent = 0;
 422
 423   AliInfo(Form(" **** analyzing event # %d  \n", fEv));
 424
 425   fTrackStore    = fMuonRecoCheck->ReconstructedTracks(fEv);
 426   if (fTrackStore->IsEmpty()) return kFALSE;
 427   fTrackRefStore = fMuonRecoCheck->ReconstructibleTracks(fEv);
 428
 429   fRunLoader->GetEvent(fEv);
 430   if (!fMFTLoader->TreeR()->GetEvent()) return kFALSE;
 431   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
 432     AliDebug(1, Form("plane %02d: nClusters = %d\n", iPlane, (fMFT->GetRecPointsList(iPlane))->GetEntries()));
 433     fMFTClusterArray[iPlane] = fMFT->GetRecPointsList(iPlane);
 434   }
 435   SeparateFrontBackClusters();
 436
 437   fRunLoader -> LoadKinematics();
 438   fStack = fRunLoader->Stack();
 439   fNextTrack = fTrackStore->CreateIterator();
 440   fMuonForwardTracks->Delete();
 441
 442   fEv++;
 443
 444   return kTRUE;
 445
 446 }
 447
 448 //======================================================================================================================================
 449
 450 Int_t AliMuonForwardTrackFinder::LoadNextTrack() {
 451
 452   fNPlanesMFTAnalyzed = 0;
 453
 454   // load next muon track from the reconstructed event
 455
 456   if (fCountRealTracksAnalyzed>fMaxNTracksToBeAnalyzed) return kFALSE;
 457
 458   if (!fCountRealTracksAnalyzed) if (!LoadNextEvent()) return kFALSE;
 459   if (fCountRealTracksAnalyzed==fMaxNTracksToBeAnalyzed) {
 460     fCountRealTracksAnalyzed++;
 461     if (!LoadNextEvent()) return kFALSE;
 462   }
 463
 464   while ( !(fMuonTrackReco = static_cast<AliMUONTrack*>(fNextTrack->Next())) ) if (!LoadNextEvent()) return kFALSE;
 465
 466   AliDebug(1, "**************************************************************************************\n");
 467   AliDebug(1, Form("***************************   MUON TRACK %3d   ***************************************\n", fCountRealTracksAnalyzedOfEvent));
 468   AliDebug(1, "**************************************************************************************\n");
 469
 470   fCountRealTracksAnalyzed++;
 471
 472   fCandidateTracks -> Delete();
 473
 474   fLabelMC = -1;
 475   fDistanceFromGoodClusterAndTrackAtLastPlane = -1.;
 476   fDistanceFromBestClusterAndTrackAtLastPlane = -1.;
 477   ResetPlanes();
 478
 479   TIter nextTrackRef(fTrackRefStore->CreateIterator());
 480   AliMUONTrack *trackRef=0;
 481
 482   // --------------------------------------- loop on MC generated tracks to find the MC reference...
 483
 484   while ( (trackRef = static_cast<AliMUONTrack*>(nextTrackRef())) ) {
 485     // number of compatible clusters between trackReco and trackRef
 486     Int_t nMatchCluster = fMuonTrackReco->FindCompatibleClusters(*trackRef, fSigmaSpectrometerCut, fIsClusterCompatible);
 487     if ( (fIsClusterCompatible[0] || fIsClusterCompatible[1] || fIsClusterCompatible[2] || fIsClusterCompatible[3]) &&   // before the dipole
 488          (fIsClusterCompatible[6] || fIsClusterCompatible[7] || fIsClusterCompatible[8] || fIsClusterCompatible[9]) &&   // after the dipole
 489          2*nMatchCluster>fMuonTrackReco->GetNClusters() ) {
 490       fMuonTrackReco->SetMCLabel(trackRef->GetUniqueID());   // MC reference has been found for trackReco!
 491       break;
 492     }
 493   }
 494
 495   // ------------------------------------- ...done!
 496
 497   fLabelMC = fMuonTrackReco->GetMCLabel();
 498
 499   Int_t motherPdg=0;
 500   if (fLabelMC>=0) {
 501     fCountRealTracksWithRefMC++;
 502     if (fStack->Particle(fLabelMC)->GetFirstMother() != -1) {
 503       motherPdg = fStack->Particle(fStack->Particle(fLabelMC)->GetFirstMother())->GetPdgCode();
 504     }
 505   }
 506
 507   CheckCurrentMuonTrackable();
 508
 509   if (fMuonTrackReco->GetMatchTrigger()) fCountRealTracksWithRefMC_andTrigger++;
 510
 511   // the track we are going to build, starting from fMuonTrackReco and adding the MFT clusters
 512   AliMuonForwardTrack *track = new ((*fCandidateTracks)[0]) AliMuonForwardTrack();
 513   track -> SetMUONTrack(fMuonTrackReco);
 514   if (fLabelMC>=0 && fStack->Particle(fLabelMC)) track->SetMCTrackRef(fStack->Particle(fLabelMC));
 515   track -> SetMCLabel(fMuonTrackReco->GetMCLabel());
 516   track -> SetMatchTrigger(fMuonTrackReco->GetMatchTrigger());
 517
 518   // track origin
 519   Double_t xVtx=-999., yVtx=-999., zVtx=-999999.;
 520   if (track->GetMCTrackRef()) {
 521     xVtx = track->GetMCTrackRef()->Vx();
 522     yVtx = track->GetMCTrackRef()->Vy();
 523     zVtx = track->GetMCTrackRef()->Vz();
 524   }
 525
 526   // track kinematics
 527   Double_t pt=-999., theta=-999., eta=-999.;
 528   if (track->GetMCTrackRef()) {
 529     pt    = track->GetMCTrackRef()->Pt();
 530     theta = track->GetMCTrackRef()->Theta();
 531     if (theta<0.) theta += TMath::Pi();
 532     eta   = track->GetMCTrackRef()->Eta();
 533   }
 534   else {
 535     AliMUONTrackParam *param = (AliMUONTrackParam*) (fMuonTrackReco->GetTrackParamAtCluster()->First());
 536     pt    = TMath::Sqrt(param->Px()*param->Px() + param->Py()*param->Py());
 537     theta = TMath::ATan(pt/param->Pz());
 538     if (theta<0.) theta += TMath::Pi();
 539     eta   = -1.*TMath::Log(TMath::Tan(0.5*theta));
 540   }
 541   // if the transverse momentum is smaller than the threshold, skip to the next track
 542   if (pt < fLowPtCut) return 3;
 543
 544   // track parameters linearly extrapolated from the first tracking station to the end of the absorber
 545   AliMUONTrackParam trackParamEndOfAbsorber(*((AliMUONTrackParam*)(fMuonTrackReco->GetTrackParamAtCluster()->First())));
 546   AliMUONTrackExtrap::ExtrapToZCov(&trackParamEndOfAbsorber, -503.);   // absorber extends from -90 to -503 cm
 547   Double_t xEndOfAbsorber = trackParamEndOfAbsorber.GetNonBendingCoor();
 548   Double_t yEndOfAbsorber = trackParamEndOfAbsorber.GetBendingCoor();
 549   Double_t rAbsorber      = TMath::Sqrt(xEndOfAbsorber*xEndOfAbsorber + yEndOfAbsorber*yEndOfAbsorber);
 550   fHistRadiusEndOfAbsorber -> Fill(rAbsorber);
 551
 552   // if the radial distance of the track at the end of the absorber is smaller than a given radius, skip to the next track
 553   if (rAbsorber < fRAbsorberCut) return 4;
 554
 555   //------------------------- NOW THE CYCLE OVER THE MFT PLANES STARTS ---------------------------------------
 556
 557   for (Int_t iPlane=fNPlanesMFT-1; iPlane>=0; iPlane--) {          // *** do not reverse the order of this cycle!!!
 558                                                                    // *** this reflects the fact that the extrapolation is performed
 559                                                                    // *** starting from the last MFT plane back to the origin
 560
 561     // --------- updating the array of tracks according to the clusters available in the i-th plane ---------
 562
 563     fNPlanesMFTAnalyzed++;
 564
 565     if (fMatchingMode==kRealMatching) {
 566       Int_t nTracksToBeAnalyzed = fCandidateTracks->GetEntriesFast();
 567       for (Int_t iTrack=0; iTrack<nTracksToBeAnalyzed; iTrack++) {
 568         fCurrentTrack = (AliMuonForwardTrack*) fCandidateTracks->UncheckedAt(iTrack);
 569         // if the old track is compatible with the new cluster, the track is updated and inserted as new track in the array
 570         // (several new tracks can be created for one old track)
 571         FindClusterInPlane(iPlane);
 572         if ((fNPlanesMFTAnalyzed-fCurrentTrack->GetNMFTClusters())>fNMaxMissingMFTClusters || fIsPlaneMandatory[iPlane]) {
 573           fCandidateTracks->Remove(fCurrentTrack);     // the old track is removed after the check;
 574         }
 575       }
 576       fCandidateTracks->Compress();
 577       if (fIsCurrentMuonTrackable) {
 578         //      fOutputQAFile->cd();
 579         fHistNTracksAfterExtrapolation[iPlane] -> Fill(fCandidateTracks->GetEntriesFast());
 580       }
 581     }
 582
 583     else if (fMatchingMode==kIdealMatching && fIsCurrentMuonTrackable) {
 584       fCurrentTrack = (AliMuonForwardTrack*) fCandidateTracks->UncheckedAt(0);
 585       AliDebug(2, Form("plane %02d: fCandidateTracks->GetEntriesFast() = %d   fCandidateTracks->UncheckedAt(0) = %p   fCurrentTrack = %p\n",
 586                        iPlane, fCandidateTracks->GetEntriesFast(), fCandidateTracks->UncheckedAt(0), fCurrentTrack));
 587       AttachGoodClusterInPlane(iPlane);
 588     }
 589
 590   }
 591
 592   // -------------------------- END OF THE CYCLE OVER THE MFT PLANES --------------------------------------------
 593
 594   AliDebug(1, "Finished cycle over planes");
 595
 596   Double_t momentum = pt * TMath::CosH(eta);
 597   fTxtTrackMomentum = new TLatex(0.10, 0.70, Form("P_{spectro} = %3.1f GeV/c", momentum));
 598
 599   if (fMatchingMode==kIdealMatching) {
 600     AliDebug(1, "Adding track to output tree...\n");
 601     fFinalBestCandidate = (AliMuonForwardTrack*) fCandidateTracks->UncheckedAt(0);
 602     AliMuonForwardTrack *newTrack = (AliMuonForwardTrack*) fCandidateTracks->UncheckedAt(0);
 603     new ((*fMuonForwardTracks)[fMuonForwardTracks->GetEntries()]) AliMuonForwardTrack(*newTrack);
 604     AliDebug(1, "...track added!\n");
 605     fCandidateTracks->Delete();
 606     fCountRealTracksAnalyzedOfEvent++;
 607     fCountRealTracksAnalyzedWithFinalCandidates++;
 608     PrintParticleHistory();
 609     FillPlanesWithTrackHistory();
 610
 611     Double_t chi2AtPlane[fNMaxPlanes] = {0};
 612     Int_t nGoodClusters = 0;
 613     Int_t nMFTClusters  = fFinalBestCandidate->GetNMFTClusters();
 614 //     Int_t nMUONClusters = fFinalBestCandidate->GetNMUONClusters();
 615     Int_t plane = 0;
 616     for (Int_t iCluster=0; iCluster<nMFTClusters; iCluster++) {
 617       while (!fFinalBestCandidate->PlaneExists(plane)) plane++;
 618       AliMFTCluster *localCluster = fFinalBestCandidate->GetMFTCluster(iCluster);
 619       chi2AtPlane[plane] = localCluster->GetLocalChi2();
 620       if (IsCorrectMatch(localCluster)) nGoodClusters++;
 621 //       Int_t nClustersGlobalTrack = nMUONClusters + (nMFTClusters-iCluster);        // Muon Spectrometer clusters + clusters in the Vertex Telescope
 622 //       Int_t ndfGlobalTrack = GetNDF(nClustersGlobalTrack);
 623 //       chi2AtPlane[plane] /= Double_t(ndfGlobalTrack);
 624       plane++;
 625     }
 626     for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
 627       fTxtTrackChi2[iPlane] = new TLatex(0.55*fRPlaneMax[fNPlanesMFT-1],
 628                                          0.90*fRPlaneMax[fNPlanesMFT-1],
 629                                          Form("#chi^{2} = %3.1f", chi2AtPlane[iPlane]));
 630     }
 631     fTxtTrackFinalChi2 = new TLatex(0.20, 0.44, Form("#chi^{2}_{final} = %3.1f", chi2AtPlane[0]));
 632
 633     if (fDrawOption) DrawPlanes();
 634     return 5;
 635   }
 636
 637   // If we have several final tracks, we must find the best candidate:
 638
 639   Int_t nFinalTracks = fCandidateTracks->GetEntriesFast();
 640   AliDebug(1, Form("nFinalTracks = %d", nFinalTracks));
 641
 642   if (nFinalTracks) fCountRealTracksAnalyzedWithFinalCandidates++;
 643
 644   Double_t theVariable_Best        = -1.;                    // variable defining the best candidate
 645   Bool_t bestCandidateExists       = kFALSE;
 646   Int_t nGoodClustersBestCandidate = 0;
 647   Int_t idBestCandidate            = 0;
 648   Double_t chi2HistoryForBestCandidate[fNMaxPlanes] = {0};  // chi2 on each plane, for the best candidate
 649   Double_t nClustersPerPlane[fNMaxPlanes] = {0};
 650   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
 651     chi2HistoryForBestCandidate[iPlane] = -1.;
 652     nClustersPerPlane[iPlane] = fMFTClusterArray[iPlane] -> GetEntries();
 653   }
 654
 655   fTxtFinalCandidates = new TLatex(0.10, 0.78, Form("N_{FinalCandidates} = %d", nFinalTracks));
 656
 657   Int_t nClustersMC = 0;
 658   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) nClustersMC += fIsGoodClusterInPlane[iPlane];
 659
 660   for (Int_t iTrack=0; iTrack<nFinalTracks; iTrack++) {
 661
 662     AliMuonForwardTrack *finalTrack = (AliMuonForwardTrack*) fCandidateTracks->UncheckedAt(iTrack);
 663
 664     Double_t chi2AtPlane[fNMaxPlanes] = {0};
 665     Int_t nGoodClusters = 0;
 666     Int_t nMFTClusters  = finalTrack->GetNMFTClusters();
 667 //     Int_t nMUONClusters = finalTrack->GetNMUONClusters();
 668
 669     Int_t plane = 0;
 670     for (Int_t iCluster=0; iCluster<nMFTClusters; iCluster++) {
 671       while (!finalTrack->PlaneExists(plane)) plane++;
 672       AliMFTCluster *localCluster = finalTrack->GetMFTCluster(iCluster);
 673       chi2AtPlane[plane] = localCluster->GetLocalChi2();
 674       if (IsCorrectMatch(localCluster)) nGoodClusters++;
 675 //       Int_t nClustersGlobalTrack = nMUONClusters + (nMFTClusters-iCluster);        // Muon Spectrometer clusters + clusters in the Vertex Telescope
 676 //       Int_t ndfGlobalTrack = GetNDF(nClustersGlobalTrack);
 677 //       chi2AtPlane[plane] /= Double_t(ndfGlobalTrack);
 678       plane++;
 679     }
 680
 681     if (fIsCurrentMuonTrackable) {
 682       //      fOutputQAFile->cd();
 683       fHistNGoodClustersForFinalTracks -> Fill(nGoodClusters);
 684     }
 685
 686     //    fOutputQAFile->cd();
 687
 688     Float_t finalCandidatesInfo[] = {Double_t(fRun),
 689                                      Double_t(fEv),
 690                                      Double_t(fCountRealTracksAnalyzedOfEvent),
 691                                      Double_t(nFinalTracks),
 692                                      Double_t(fLabelMC>=0),
 693                                      xVtx, yVtx, zVtx,
 694                                      motherPdg,
 695                                      Double_t(fMuonTrackReco->GetMatchTrigger()),
 696                                      Double_t(nClustersMC),
 697                                      Double_t(nGoodClusters),
 698                                      pt, theta, eta,
 699                                      chi2AtPlane[0],
 700                                      chi2AtPlane[1],
 701                                      chi2AtPlane[2],
 702                                      chi2AtPlane[3],
 703                                      chi2AtPlane[4],
 704                                      chi2AtPlane[5],
 705                                      chi2AtPlane[6],
 706                                      chi2AtPlane[7],
 707                                      chi2AtPlane[8]};
 708
 709     fNtuFinalCandidates -> Fill(finalCandidatesInfo);
 710
 711     // now comparing the tracks with various criteria, in order to find the best one
 712
 713     Double_t theVariable = 0.;
 714 //     theVariable = chi2AtPlane[0];
 715     for (Int_t iCluster=0; iCluster<nMFTClusters; iCluster++) theVariable += chi2AtPlane[iCluster];
 716     theVariable /= Double_t(nMFTClusters);
 717
 718
 719     if (theVariable<theVariable_Best || theVariable_Best<0.) {
 720       nGoodClustersBestCandidate = nGoodClusters;
 721       for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) chi2HistoryForBestCandidate[iPlane] = chi2AtPlane[iPlane];
 722       theVariable_Best = theVariable;
 723       fTxtTrackFinalChi2 = new TLatex(0.20, 0.44, Form("#chi^{2}_{final} = %3.1f", chi2HistoryForBestCandidate[0]));
 724       for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
 725         fTxtTrackChi2[iPlane] = new TLatex(0.55*fRPlaneMax[fNPlanesMFT-1],
 726                                            0.90*fRPlaneMax[fNPlanesMFT-1],
 727                                            Form("#chi^{2} = %3.1f", chi2AtPlane[iPlane]));
 728       }
 729       idBestCandidate = iTrack;
 730       bestCandidateExists=kTRUE;
 731     }
 732
 733     // ----------------------------------------------------------
 734
 735   }
 736
 737   if (nFinalTracks) {
 738     fFinalBestCandidate = (AliMuonForwardTrack*) fCandidateTracks->UncheckedAt(idBestCandidate);
 739     AliInfo(Form("fFinalBestCandidate->GetNMFTClusters() = %d\n",  fFinalBestCandidate->GetNMFTClusters()));
 740     PrintParticleHistory();
 741     FillPlanesWithTrackHistory();
 742     AliMuonForwardTrack *newTrack = (AliMuonForwardTrack*) fCandidateTracks->UncheckedAt(idBestCandidate);
 743     newTrack -> SetNWrongClustersMC(newTrack->GetNMFTClusters() - nGoodClustersBestCandidate);
 744     newTrack -> SetTrackMCId(fRun*100000+fEv*1000+fCountRealTracksAnalyzedOfEvent);
 745     new ((*fMuonForwardTracks)[fMuonForwardTracks->GetEntries()]) AliMuonForwardTrack(*newTrack);
 746   }
 747
 748   //  fOutputQAFile->cd();
 749
 750   Float_t finalBestCandidatesInfo[] = {Double_t(fRun),
 751                                        Double_t(fEv),
 752                                        Double_t(fCountRealTracksAnalyzedOfEvent),
 753                                        Double_t(nFinalTracks),
 754                                        Double_t(fLabelMC>=0),
 755                                        xVtx, yVtx, zVtx,
 756                                        motherPdg,
 757                                        Double_t(fMuonTrackReco->GetMatchTrigger()),
 758                                        Double_t(nClustersMC),
 759                                        Double_t(nGoodClustersBestCandidate),
 760                                        pt, theta, eta,
 761                                        chi2HistoryForBestCandidate[0],
 762                                        chi2HistoryForBestCandidate[1],
 763                                        chi2HistoryForBestCandidate[2],
 764                                        chi2HistoryForBestCandidate[3],
 765                                        chi2HistoryForBestCandidate[4],
 766                                        chi2HistoryForBestCandidate[5],
 767                                        chi2HistoryForBestCandidate[6],
 768                                        chi2HistoryForBestCandidate[7],
 769                                        chi2HistoryForBestCandidate[8],
 770                                        nClustersPerPlane[0],
 771                                        nClustersPerPlane[1],
 772                                        nClustersPerPlane[2],
 773                                        nClustersPerPlane[3],
 774                                        nClustersPerPlane[4],
 775                                        nClustersPerPlane[5],
 776                                        nClustersPerPlane[6],
 777                                        nClustersPerPlane[7],
 778                                        nClustersPerPlane[8]};
 779
 780   fNtuFinalBestCandidates -> Fill(finalBestCandidatesInfo);
 781
 782   if (fDrawOption && bestCandidateExists) {
 783     fTxtTrackGoodClusters = new TLatex(0.20, 0.51, Form("N_{GoodClusters} = %d", nGoodClustersBestCandidate));
 784     DrawPlanes();
 785   }
 786
 787   // -------------------------------------------------------------------------------------------
 788
 789   fCandidateTracks->Delete();
 790   fFinalBestCandidate = NULL;
 791
 792   fCountRealTracksAnalyzedOfEvent++;
 793
 794   return 5;
 795
 796 }
 797
 798 //===========================================================================================================================================
 799
 800 void AliMuonForwardTrackFinder::FindClusterInPlane(Int_t planeId) {
 801
 802   AliDebug(2, Form(">>>> executing AliMuonForwardTrackFinder::FindClusterInPlane(%d)\n", planeId));
 803
 804   // !!!!!!!!! coordinates and errors on the interaction vertex should be taken from the event itself (ITS) if available
 805
 806   // propagate track to plane #planeId (both to front and back active sensors)
 807   // look for compatible clusters
 808   // update TrackParam at found cluster (if any) using Kalman Filter
 809
 810   AliMUONTrackParam currentParamFront, currentParamBack, currentParamForResearchFront, currentParamForResearchBack;
 811
 812   if (planeId == fNPlanesMFT-1) {      // last plane of the telecope
 813     currentParamFront = (*((AliMUONTrackParam*)(fMuonTrackReco->GetTrackParamAtCluster()->First())));
 814     currentParamBack  = (*((AliMUONTrackParam*)(fMuonTrackReco->GetTrackParamAtCluster()->First())));
 815     currentParamForResearchFront = currentParamFront;
 816     currentParamForResearchBack  = currentParamBack;
 817     Double_t xExtrap = gRandom->Gaus(0,fVertexErrorX);
 818     Double_t yExtrap = gRandom->Gaus(0,fVertexErrorY);
 819     Double_t zExtrap = gRandom->Gaus(0,fVertexErrorZ);
 820     if (fBransonCorrection) {
 821       AliMUONTrackExtrap::ExtrapToVertex(&currentParamFront, xExtrap, yExtrap, zExtrap, fVertexErrorX, fVertexErrorY);
 822       AliMUONTrackExtrap::ExtrapToVertex(&currentParamBack,  xExtrap, yExtrap, zExtrap, fVertexErrorX, fVertexErrorY);
 823     }
 824     else {
 825       AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(&currentParamFront, zExtrap);
 826       AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(&currentParamBack,  zExtrap);
 827     }
 828     AliMUONTrackExtrap::ExtrapToVertex(&currentParamForResearchFront, xExtrap, yExtrap, zExtrap, fVertexErrorX, fVertexErrorY);
 829     AliMUONTrackExtrap::ExtrapToVertex(&currentParamForResearchBack,  xExtrap, yExtrap, zExtrap, fVertexErrorX, fVertexErrorY);
 830   }
 831   else {          // MFT planes others than the last one: mult. scattering correction because of the upstream MFT planes is performed
 832     currentParamFront = (*((AliMUONTrackParam*)(fCurrentTrack->GetTrackParamAtCluster()->First())));
 833     currentParamBack  = (*((AliMUONTrackParam*)(fCurrentTrack->GetTrackParamAtCluster()->First())));
 834     currentParamForResearchFront = currentParamFront;
 835     currentParamForResearchBack  = currentParamBack;
 836     AliMUONTrackExtrap::AddMCSEffect(&currentParamFront,           (fSegmentation->GetPlane(planeId+1)->GetEquivalentSilicon()+
 837                                                                     fSegmentation->GetPlane(planeId)->GetEquivalentSiliconBeforeFront())/fRadLengthSi,-1.);
 838     AliMUONTrackExtrap::AddMCSEffect(&currentParamForResearchFront,(fSegmentation->GetPlane(planeId+1)->GetEquivalentSilicon()+
 839                                                                     fSegmentation->GetPlane(planeId)->GetEquivalentSiliconBeforeFront())/fRadLengthSi,-1.);
 840     AliMUONTrackExtrap::AddMCSEffect(&currentParamBack,            (fSegmentation->GetPlane(planeId+1)->GetEquivalentSilicon()+
 841                                                                     fSegmentation->GetPlane(planeId)->GetEquivalentSiliconBeforeBack())/fRadLengthSi,-1.);
 842     AliMUONTrackExtrap::AddMCSEffect(&currentParamForResearchBack, (fSegmentation->GetPlane(planeId+1)->GetEquivalentSilicon()+
 843                                                                     fSegmentation->GetPlane(planeId)->GetEquivalentSiliconBeforeBack())/fRadLengthSi,-1.);
 844   }
 845   // for all planes: extrapolation to the Z of the plane
 846   AliMUONTrackExtrap::ExtrapToZCov(&currentParamFront,            -1.*fSegmentation->GetPlane(planeId)->GetZCenterActiveFront());
 847   AliMUONTrackExtrap::ExtrapToZCov(&currentParamForResearchFront, -1.*fSegmentation->GetPlane(planeId)->GetZCenterActiveFront());
 848   AliMUONTrackExtrap::ExtrapToZCov(&currentParamBack,             -1.*fSegmentation->GetPlane(planeId)->GetZCenterActiveBack());
 849   AliMUONTrackExtrap::ExtrapToZCov(&currentParamForResearchBack,  -1.*fSegmentation->GetPlane(planeId)->GetZCenterActiveBack());
 850
 851   //---------------------------------------------------------------------------------------
 852
 853   TMatrixD covFront(5,5); covFront = currentParamForResearchFront.GetCovariances();
 854   TMatrixD covBack(5,5);  covBack  = currentParamForResearchBack.GetCovariances();
 855
 856   Double_t squaredError_X_Front = covFront(0,0);
 857   Double_t squaredError_Y_Front = covFront(2,2);
 858   Double_t squaredError_X_Back  = covBack(0,0);
 859   Double_t squaredError_Y_Back  = covBack(2,2);
 860
 861   Double_t corrFact = 1.0;
 862
 863   Double_t researchRadiusFront = TMath::Sqrt(squaredError_X_Front + squaredError_Y_Front);
 864   Double_t researchRadiusBack  = TMath::Sqrt(squaredError_X_Back  + squaredError_Y_Back);
 865   if (0.5*(researchRadiusFront+researchRadiusBack)<fMinResearchRadiusAtPlane[planeId]) {
 866     corrFact = fMinResearchRadiusAtPlane[planeId]/(0.5*(researchRadiusFront+researchRadiusBack));
 867   }
 868   if (fIsCurrentMuonTrackable) {
 869     //    fOutputQAFile->cd();
 870     fHistResearchRadius[planeId] -> Fill(0.5*(researchRadiusFront+researchRadiusBack));
 871   }
 872
 873   Double_t position_X_Front = currentParamForResearchFront.GetNonBendingCoor();
 874   Double_t position_Y_Front = currentParamForResearchFront.GetBendingCoor();
 875   Double_t position_X_Back = currentParamForResearchBack.GetNonBendingCoor();
 876   Double_t position_Y_Back = currentParamForResearchBack.GetBendingCoor();
 877   Double_t radialPositionOfTrackFront = TMath::Sqrt(position_X_Front*position_X_Front + position_Y_Front*position_Y_Front);
 878   Double_t radialPositionOfTrackBack  = TMath::Sqrt(position_X_Back*position_X_Back   + position_Y_Back*position_Y_Back);
 879
 880   //---------------------------------------------------------------------------------------
 881
 882   Double_t chi2cut = 2.*fSigmaClusterCut*fSigmaClusterCut;     // depends on the number of variables (here, 2)
 883
 884   // Analyizing the clusters: FRONT ACTIVE ELEMENTS
 885
 886   Int_t nClustersFront = fMFTClusterArrayFront[planeId]->GetEntries();
 887   AliDebug(2, Form("There are %3d clusters in plane %02d FRONT\n", nClustersFront, planeId));
 888
 889   for (Int_t iCluster=0; iCluster<nClustersFront; iCluster++) {
 890
 891     Bool_t isGoodChi2 = kFALSE;
 892
 893     AliMFTCluster *cluster = (AliMFTCluster*) fMFTClusterArrayFront[planeId]->At(iCluster);
 894     Double_t chi2 = (1./(corrFact*corrFact)) * TryOneCluster(currentParamForResearchFront, cluster);     // describes the compatibility between the track and the cluster
 895     if (chi2<chi2cut) isGoodChi2 = kTRUE;
 896
 897     Double_t radialPositionOfClusterFront = TMath::Sqrt(cluster->GetX()*cluster->GetX() + cluster->GetY()*cluster->GetY());
 898     if (planeId == fNPlanesMFT-1) {
 899       if (TMath::Abs(radialPositionOfTrackFront-radialPositionOfClusterFront)<fDistanceFromBestClusterAndTrackAtLastPlane ||
 900           fDistanceFromBestClusterAndTrackAtLastPlane<0.) {
 901         fDistanceFromBestClusterAndTrackAtLastPlane = TMath::Abs(radialPositionOfTrackFront-radialPositionOfClusterFront);
 902       }
 903       if (IsCorrectMatch(cluster)) {
 904         fDistanceFromGoodClusterAndTrackAtLastPlane = TMath::Abs(radialPositionOfTrackFront-radialPositionOfClusterFront);
 905       }
 906     }
 907
 908     if (fIsCurrentMuonTrackable) {
 909       //      fOutputQAFile->cd();
 910       if (IsCorrectMatch(cluster)) fHistChi2Cluster_GoodCluster[planeId]->Fill(chi2/2.);     //  chi2/ndf
 911       else                         fHistChi2Cluster_BadCluster[planeId] ->Fill(chi2/2.);     //  chi2/ndf
 912     }
 913
 914     if (isGoodChi2) {
 915       AliDebug(3, Form("accepting cluster: chi2=%f (cut = %f)\n", chi2, chi2cut));
 916       AliMuonForwardTrack *newTrack = new ((*fCandidateTracks)[fCandidateTracks->GetEntriesFast()]) AliMuonForwardTrack(*fCurrentTrack);
 917       newTrack->AddTrackParamAtMFTCluster(currentParamFront, *cluster);    // creating new track param and attaching the cluster
 918       AliDebug(1, Form("After plane %02d: newTrack->GetNMFTClusters() = %d (fCurrentTrack->GetNMFTClusters() = %d)",
 919                        planeId, newTrack->GetNMFTClusters(), fCurrentTrack->GetNMFTClusters()));
 920       newTrack->SetPlaneExists(planeId);
 921       AliDebug(2, Form("current muon is trackable: %d\n", fIsCurrentMuonTrackable));
 922       if (fIsCurrentMuonTrackable) {
 923         Double_t newGlobalChi2 = ((AliMUONTrackParam*) newTrack->GetTrackParamAtCluster()->First())->GetTrackChi2();
 924         AliDebug(2, Form("new chi2 = %f (= %f)\n", newGlobalChi2, newTrack->GetMFTCluster(0)->GetTrackChi2()));
 925         Int_t nClustersGlobalTrack = newTrack->GetNMUONClusters() + newTrack->GetNMFTClusters();        // Muon Spectrometer clusters + clusters in the Vertex Telescope
 926         Int_t ndfGlobalTrack = GetNDF(nClustersGlobalTrack);
 927         //      fOutputQAFile->cd();
 928         if (IsCorrectMatch(cluster)) fHistGlobalChi2AtPlaneFor_GOOD_CandidatesOfTrackableMuons[planeId]->Fill(newGlobalChi2/Double_t(ndfGlobalTrack));
 929         else                         fHistGlobalChi2AtPlaneFor_BAD_CandidatesOfTrackableMuons[planeId] ->Fill(newGlobalChi2/Double_t(ndfGlobalTrack));
 930       }
 931       fGrMFTPlane[kClustersGoodChi2][planeId] -> SetPoint(fGrMFTPlane[kClustersGoodChi2][planeId]->GetN(), cluster->GetX(), cluster->GetY());
 932     }
 933     else AliDebug(3, Form("discarding cluster: chi2=%f (cut = %f)\n", chi2, chi2cut));
 934
 935   }
 936
 937   // Analyizing the clusters: BACK ACTIVE ELEMENTS
 938
 939   Int_t nClustersBack = fMFTClusterArrayBack[planeId]->GetEntries();
 940   AliDebug(2, Form("There are %3d clusters in plane %02d BACK\n", nClustersBack, planeId));
 941
 942   for (Int_t iCluster=0; iCluster<nClustersBack; iCluster++) {
 943
 944     Bool_t isGoodChi2 = kFALSE;
 945
 946     AliMFTCluster *cluster = (AliMFTCluster*) fMFTClusterArrayBack[planeId]->At(iCluster);
 947     Double_t chi2 = (1./(corrFact*corrFact)) * TryOneCluster(currentParamForResearchBack, cluster);     // describes the compatibility between the track and the cluster
 948     if (chi2<chi2cut) isGoodChi2 = kTRUE;
 949
 950     Double_t radialPositionOfClusterBack = TMath::Sqrt(cluster->GetX()*cluster->GetX() + cluster->GetY()*cluster->GetY());
 951     if (planeId == fNPlanesMFT-1) {
 952       if (TMath::Abs(radialPositionOfTrackBack-radialPositionOfClusterBack)<fDistanceFromBestClusterAndTrackAtLastPlane ||
 953           fDistanceFromBestClusterAndTrackAtLastPlane<0.) {
 954         fDistanceFromBestClusterAndTrackAtLastPlane = TMath::Abs(radialPositionOfTrackBack-radialPositionOfClusterBack);
 955       }
 956       if (IsCorrectMatch(cluster)) {
 957         fDistanceFromGoodClusterAndTrackAtLastPlane = TMath::Abs(radialPositionOfTrackBack-radialPositionOfClusterBack);
 958       }
 959     }
 960
 961     if (fIsCurrentMuonTrackable) {
 962       //      fOutputQAFile->cd();
 963       if (IsCorrectMatch(cluster)) fHistChi2Cluster_GoodCluster[planeId]->Fill(chi2/2.);     //  chi2/ndf
 964       else                         fHistChi2Cluster_BadCluster[planeId] ->Fill(chi2/2.);     //  chi2/ndf
 965     }
 966
 967     if (isGoodChi2) {
 968       AliDebug(3,Form("accepting cluster: chi2=%f (cut = %f)\n", chi2, chi2cut));
 969       AliMuonForwardTrack *newTrack = new ((*fCandidateTracks)[fCandidateTracks->GetEntriesFast()]) AliMuonForwardTrack(*fCurrentTrack);
 970       newTrack->AddTrackParamAtMFTCluster(currentParamBack, *cluster);    // creating new track param and attaching the cluster
 971       AliDebug(1, Form("After plane %02d: newTrack->GetNMFTClusters() = %d (fCurrentTrack->GetNMFTClusters() = %d)",
 972                        planeId, newTrack->GetNMFTClusters(), fCurrentTrack->GetNMFTClusters()));
 973       newTrack->SetPlaneExists(planeId);
 974       AliDebug(2, Form("current muon is trackable: %d\n", fIsCurrentMuonTrackable));
 975       if (fIsCurrentMuonTrackable) {
 976         Double_t newGlobalChi2 = ((AliMUONTrackParam*) newTrack->GetTrackParamAtCluster()->First())->GetTrackChi2();
 977         AliDebug(2, Form("new chi2 = %f (= %f)\n", newGlobalChi2, newTrack->GetMFTCluster(0)->GetTrackChi2()));
 978         Int_t nClustersGlobalTrack = newTrack->GetNMUONClusters() + newTrack->GetNMFTClusters();        // Muon Spectrometer clusters + clusters in the Vertex Telescope
 979         Int_t ndfGlobalTrack = GetNDF(nClustersGlobalTrack);
 980         //      fOutputQAFile->cd();
 981         if (IsCorrectMatch(cluster)) fHistGlobalChi2AtPlaneFor_GOOD_CandidatesOfTrackableMuons[planeId]->Fill(newGlobalChi2/Double_t(ndfGlobalTrack));
 982         else                         fHistGlobalChi2AtPlaneFor_BAD_CandidatesOfTrackableMuons[planeId] ->Fill(newGlobalChi2/Double_t(ndfGlobalTrack));
 983       }
 984       fGrMFTPlane[kClustersGoodChi2][planeId] -> SetPoint(fGrMFTPlane[kClustersGoodChi2][planeId]->GetN(), cluster->GetX(), cluster->GetY());
 985     }
 986     else AliDebug(3,Form("discarding cluster: chi2=%f (cut = %f)\n", chi2, chi2cut));
 987
 988   }
 989
 990   //---------------------------------------------------------------------------------------------
 991
 992   if (planeId == fNPlanesMFT-1) {
 993     if (fIsCurrentMuonTrackable && fDistanceFromGoodClusterAndTrackAtLastPlane>0.) {
 994       //      fOutputQAFile->cd();
 995       fHistDistanceGoodClusterFromTrackMinusDistanceBestClusterFromTrackAtLastPlane -> Fill(TMath::Abs(fDistanceFromBestClusterAndTrackAtLastPlane-
 996                                                                                                        fDistanceFromGoodClusterAndTrackAtLastPlane));
 997       fHistDistanceGoodClusterFromTrackAtLastPlane -> Fill(fDistanceFromGoodClusterAndTrackAtLastPlane);
 998     }
 999   }
1000
1001 }
1002
1003 //==========================================================================================================================================
1004
1005 void AliMuonForwardTrackFinder::AttachGoodClusterInPlane(Int_t planeId) {
1006
1007   AliDebug(1, Form(">>>> executing AliMuonForwardTrackFinder::AttachGoodClusterInPlane(%d)\n", planeId));
1008
1009   AliMUONTrackParam currentParamFront, currentParamBack;
1010
1011   if (planeId == fNPlanesMFT-1) {      // last plane of the telecope
1012     currentParamFront = (*((AliMUONTrackParam*)(fMuonTrackReco->GetTrackParamAtCluster()->First())));
1013     currentParamBack  = (*((AliMUONTrackParam*)(fMuonTrackReco->GetTrackParamAtCluster()->First())));
1014     AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(&currentParamFront, 0.);
1015     AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(&currentParamBack,  0.);
1016   }
1017   else {          // MFT planes others than the last one: mult. scattering correction because of the upstream MFT planes is performed
1018     AliDebug(2, Form("fCurrentTrack = %p\n", fCurrentTrack));
1019     currentParamFront = (*((AliMUONTrackParam*)(fCurrentTrack->GetTrackParamAtCluster()->First())));
1020     currentParamBack  = (*((AliMUONTrackParam*)(fCurrentTrack->GetTrackParamAtCluster()->First())));
1021     AliMUONTrackExtrap::AddMCSEffect(&currentParamFront, (fSegmentation->GetPlane(planeId+1)->GetEquivalentSilicon()+
1022                                                           fSegmentation->GetPlane(planeId)->GetEquivalentSiliconBeforeFront())/fRadLengthSi,-1.);
1023     AliMUONTrackExtrap::AddMCSEffect(&currentParamBack,  (fSegmentation->GetPlane(planeId+1)->GetEquivalentSilicon()+
1024                                                           fSegmentation->GetPlane(planeId)->GetEquivalentSiliconBeforeBack())/fRadLengthSi,-1.);
1025   }
1026   // for all planes: linear extrapolation to the Z of the plane
1027   AliMUONTrackExtrap::ExtrapToZCov(&currentParamFront, -1.*fSegmentation->GetPlane(planeId)->GetZCenterActiveFront());
1028   AliMUONTrackExtrap::ExtrapToZCov(&currentParamBack,  -1.*fSegmentation->GetPlane(planeId)->GetZCenterActiveBack());
1029
1030   Bool_t goodClusterFound = kFALSE;
1031
1032   // Analyizing the clusters: FRONT ACTIVE ELEMENTS
1033
1034   Int_t nClustersFront = fMFTClusterArrayFront[planeId]->GetEntries();
1035
1036   AliDebug(1, Form("nClustersFront = %d\n", nClustersFront));
1037   for (Int_t iCluster=0; iCluster<nClustersFront; iCluster++) {
1038     AliMFTCluster *cluster = (AliMFTCluster*) fMFTClusterArrayFront[planeId]->UncheckedAt(iCluster);
1039     AliDebug(2, Form("checking cluster %02d of %02d: cluter=%p, fCurrentTrack=%p\n", iCluster, nClustersFront, cluster, fCurrentTrack));
1040     if (IsCorrectMatch(cluster)) {
1041       fCurrentTrack->AddTrackParamAtMFTCluster(currentParamFront, *cluster);  // creating new track param and attaching the cluster
1042       fCurrentTrack->SetPlaneExists(planeId);
1043       goodClusterFound = kTRUE;
1044       break;
1045     }
1046   }
1047
1048   if (goodClusterFound) return;
1049
1050   // Analyizing the clusters: BACK ACTIVE ELEMENTS
1051
1052   Int_t nClustersBack = fMFTClusterArrayBack[planeId]->GetEntries();
1053
1054   AliDebug(1, Form("nClustersBack = %d\n", nClustersBack));
1055   for (Int_t iCluster=0; iCluster<nClustersBack; iCluster++) {
1056     AliMFTCluster *cluster = (AliMFTCluster*) fMFTClusterArrayBack[planeId]->UncheckedAt(iCluster);
1057     AliDebug(2,Form("checking cluster %02d of %02d: cluter=%p, fCurrentTrack=%p\n", iCluster, nClustersBack, cluster, fCurrentTrack));
1058     if (IsCorrectMatch(cluster)) {
1059       fCurrentTrack->AddTrackParamAtMFTCluster(currentParamBack, *cluster);  // creating new track param and attaching the cluster
1060       fCurrentTrack->SetPlaneExists(planeId);
1061       goodClusterFound = kTRUE;
1062       break;
1063     }
1064   }
1065
1066 }
1067
1068 //==========================================================================================================================================
1069
1070 void AliMuonForwardTrackFinder::CheckCurrentMuonTrackable() {
1071
1072   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
1073     fIsGoodClusterInPlane[iPlane] = kFALSE;
1074     Int_t nClusters = fMFTClusterArray[iPlane]->GetEntriesFast();
1075     for (Int_t iCluster=0; iCluster<nClusters; iCluster++) {
1076       AliMFTCluster *cluster = (AliMFTCluster*) fMFTClusterArray[iPlane]->At(iCluster);
1077       for (Int_t iTrack=0; iTrack<cluster->GetNMCTracks(); iTrack++) {
1078         if (cluster->GetMCLabel(iTrack)==fLabelMC) {
1079           fIsGoodClusterInPlane[iPlane] = kTRUE;
1080           break;
1081         }
1082       }
1083     }
1084   }
1085
1086   fIsCurrentMuonTrackable = kTRUE;
1087   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) fIsCurrentMuonTrackable = (fIsCurrentMuonTrackable&&fIsGoodClusterInPlane[iPlane]);
1088
1089 }
1090
1091 //==========================================================================================================================================
1092
1093 void AliMuonForwardTrackFinder::FillPlanesWithTrackHistory() {
1094
1095   // Fill planes with the clusters
1096
1097   Int_t cluster = 0;
1098   AliDebug(2, Form("fFinalBestCandidate->GetNMFTClusters() = %d\n",  fFinalBestCandidate->GetNMFTClusters()));
1099   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
1100     if (fFinalBestCandidate->PlaneExists(iPlane)) {
1101       AliMFTCluster *trackCluster = fFinalBestCandidate->GetMFTCluster(cluster++);
1102       fGrMFTPlane[kClusterOfTrack][iPlane] -> SetPoint(fGrMFTPlane[kClusterOfTrack][iPlane]->GetN(), trackCluster->GetX(), trackCluster->GetY());
1103     }
1104     Int_t nClusters = fMFTClusterArray[iPlane]->GetEntriesFast();
1105     for (Int_t iCluster=0; iCluster<nClusters; iCluster++) {
1106       AliMFTCluster *myCluster = (AliMFTCluster*) fMFTClusterArray[iPlane]->UncheckedAt(iCluster);
1107       fGrMFTPlane[kAllClusters][iPlane] -> SetPoint(fGrMFTPlane[kAllClusters][iPlane]->GetN(), myCluster->GetX(), myCluster->GetY());
1108       if (IsCorrectMatch(myCluster)) {
1109         fGrMFTPlane[kClusterCorrectMC][iPlane] -> SetPoint(fGrMFTPlane[kClusterCorrectMC][iPlane]->GetN(), myCluster->GetX(), myCluster->GetY());
1110       }
1111     }
1112   }
1113
1114 }
1115
1116 //======================================================================================================================================
1117
1118 Bool_t AliMuonForwardTrackFinder::IsCorrectMatch(AliMFTCluster *cluster) {
1119
1120   Bool_t result = kFALSE;
1121
1122   // check if the cluster belongs to the correct MC track
1123
1124   for (Int_t iTrack=0; iTrack<cluster->GetNMCTracks(); iTrack++) {
1125     if (cluster->GetMCLabel(iTrack)==fLabelMC) {
1126       result = kTRUE;
1127       break;
1128     }
1129   }
1130
1131   AliDebug(2,Form("returning %d\n", result));
1132
1133   return result;
1134
1135 }
1136
1137 //======================================================================================================================================
1138
1139 Double_t AliMuonForwardTrackFinder::TryOneCluster(const AliMUONTrackParam &trackParam, AliMFTCluster *cluster) {
1140
1141   // Test the compatibility between the track and the cluster (using trackParam's covariance matrix):
1142   // return the corresponding Chi2
1143   // assume the track parameters are given at the Z of the cluster
1144
1145   // Set differences between trackParam and cluster in the bending and non bending directions
1146   Double_t dX = cluster->GetX() - trackParam.GetNonBendingCoor();
1147   Double_t dY = cluster->GetY() - trackParam.GetBendingCoor();
1148   AliDebug(3,Form("dX = %f, dY = %f\n", dX, dY));
1149
1150   // Calculate errors and covariances
1151   const TMatrixD& kParamCov = trackParam.GetCovariances();
1152   Double_t sigmaX2 = kParamCov(0,0) + cluster->GetErrX2();
1153   Double_t sigmaY2 = kParamCov(2,2) + cluster->GetErrY2();
1154   AliDebug(3, Form("dX2 = %f, dY2 = %f\n", sigmaX2, sigmaY2));
1155   Double_t covXY   = kParamCov(0,2);
1156   Double_t det     = sigmaX2 * sigmaY2 - covXY * covXY;
1157
1158   // Compute chi2
1159   if (det==0.) return 1.e10;
1160   return (dX*dX*sigmaY2 + dY*dY*sigmaX2 - 2.*dX*dY*covXY) / det;
1161
1162 }
1163
1164 //=========================================================================================================================================
1165
1166 void AliMuonForwardTrackFinder::SeparateFrontBackClusters() {
1167
1168   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
1169     fMFTClusterArrayFront[iPlane]->Delete();
1170     fMFTClusterArrayBack[iPlane] ->Delete();
1171     for (Int_t iCluster=0; iCluster<fMFTClusterArray[iPlane]->GetEntries(); iCluster++) {
1172       AliMFTCluster *cluster = (AliMFTCluster*) fMFTClusterArray[iPlane]->At(iCluster);
1173       if (TMath::Abs(cluster->GetZ())<TMath::Abs(fSegmentation->GetPlane(iPlane)->GetZCenter())) {
1174         new ((*fMFTClusterArrayFront[iPlane])[fMFTClusterArrayFront[iPlane]->GetEntries()]) AliMFTCluster(*cluster);
1175       }
1176       else {
1177         new ((*fMFTClusterArrayBack[iPlane])[fMFTClusterArrayBack[iPlane]->GetEntries()]) AliMFTCluster(*cluster);
1178       }
1179     }
1180   }
1181
1182 }
1183
1184 //=========================================================================================================================================
1185
1186 Int_t AliMuonForwardTrackFinder::GetNDF(Int_t nClusters) {
1187
1188   // the same definition as in AliMUONTrack is implemented, since here we just add more clusters to the Muon track
1189
1190   Int_t ndf = 2 * nClusters - 5;
1191   return (ndf > 0) ? ndf : 0;
1192
1193 }
1194
1195 //============================================================================================================================================
1196
1197 void AliMuonForwardTrackFinder::BookHistos() {
1198
1199   const Int_t nMaxNewTracks[]  = {150,     200,   250, 600, 1000};
1200   const Double_t radiusPlane[] = {0.010, 0.010, 0.050, 0.5,  1.5};
1201
1202   fHistRadiusEndOfAbsorber = new TH1D("hRadiusEndOfAbsorber", "Track radial distance at the end of the absorber",  1000, 0, 100.);
1203
1204   fHistNGoodClustersForFinalTracks = new TH1D("hNGoodClustersForFinalTracks", "Number of Good Clusters per Final Track", 20, -0.25, 9.75);
1205
1206   fHistDistanceGoodClusterFromTrackAtLastPlane = new TH1D("hDistanceGoodClusterFromTrackAtLastPlane",
1207                                                           "Distance of MC Good Cluster from Track in last MFT plane", 200, 0., 2.);
1208
1209   fHistDistanceGoodClusterFromTrackMinusDistanceBestClusterFromTrackAtLastPlane =
1210     new TH1D("hDistanceGoodClusterFromTrackMinusDistanceBestClusterFromTrackAtLastPlane",
1211              "Good Cluster distance from track - Best Cluster distance from track in last MFT plane", 200, 0., 2.);
1212
1213   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
1214
1215     fHistNTracksAfterExtrapolation[iPlane] = new TH1D(Form("hNTracksAfterExtrapolation_pl%02d", iPlane),
1216                                                       Form("Number of Candidates after analysis of MFT plane %02d", iPlane),
1217                                                       nMaxNewTracks[iPlane], -0.5, nMaxNewTracks[iPlane]-0.5);
1218
1219     fHistResearchRadius[iPlane] = new TH1D(Form("hResearchRadius_pl%02d", iPlane),
1220                                            Form("Research Radius for candidate clusters in MFT plane %02d", iPlane),
1221                                            1000, 0., radiusPlane[iPlane]);
1222
1223     fHistChi2Cluster_GoodCluster[iPlane] = new TH1D(Form("hChi2Cluster_GoodCluster_pl%02d", iPlane),
1224                                                     Form("#chi^{2}_{clust} for Good clusters in MFT plane %02d", iPlane),
1225                                                     100, 0., 15.);
1226
1227     fHistChi2Cluster_BadCluster[iPlane] = new TH1D(Form("hChi2Cluster_BadCluster_pl%02d", iPlane),
1228                                                    Form("#chi^{2}_{clust} for Bad clusters in MFT plane %02d", iPlane),
1229                                                    100, 0., 15.);
1230
1231     fHistGlobalChi2AtPlaneFor_GOOD_CandidatesOfTrackableMuons[iPlane] = new TH1D(Form("fHistGlobalChi2AtPlaneFor_GOOD_CandidatesOfTrackableMuons_pl%02d", iPlane),
1232                                                                                  Form("#chi^{2}/ndf at plane %d for GOOD candidates of trackable muons",iPlane),
1233                                                                                  100, 0., 15.);
1234
1235     fHistGlobalChi2AtPlaneFor_BAD_CandidatesOfTrackableMuons[iPlane] = new TH1D(Form("fHistGlobalChi2AtPlaneFor_BAD_CandidatesOfTrackableMuons_pl%02d", iPlane),
1236                                                                                 Form("#chi^{2}/ndf at plane %d for BAD candidates of trackable muons",iPlane),
1237                                                                                 100, 0., 15.);
1238
1239   }
1240
1241   //------------------------------------------
1242
1243   fHistRadiusEndOfAbsorber          -> Sumw2();
1244   fHistNGoodClustersForFinalTracks  -> Sumw2();
1245
1246   fHistDistanceGoodClusterFromTrackAtLastPlane                                  -> Sumw2();
1247   fHistDistanceGoodClusterFromTrackMinusDistanceBestClusterFromTrackAtLastPlane -> Sumw2();
1248
1249   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
1250
1251     fHistNTracksAfterExtrapolation[iPlane] -> Sumw2();
1252     fHistResearchRadius[iPlane]            -> Sumw2();
1253
1254     fHistChi2Cluster_GoodCluster[iPlane]        -> Sumw2();
1255     fHistChi2Cluster_BadCluster[iPlane]         -> Sumw2();
1256
1257     fHistGlobalChi2AtPlaneFor_GOOD_CandidatesOfTrackableMuons[iPlane] -> Sumw2();
1258     fHistGlobalChi2AtPlaneFor_BAD_CandidatesOfTrackableMuons[iPlane]  -> Sumw2();
1259
1260   }
1261
1262   fNtuFinalCandidates     = new TNtuple("ntuFinalCandidates",     "Final Candidates (ALL)", "run:event:muonTrack:nFinalCandidates:MCTrackRefExists:xVtx:yVtx:zVtx:motherPdg:triggerMatch:nClustersMC:nGoodClusters:pt:theta:eta:chi2AtPlane0:chi2AtPlane1:chi2AtPlane2:chi2AtPlane3:chi2AtPlane4:chi2AtPlane5:chi2AtPlane6:chi2AtPlane7:chi2AtPlane8");
1263
1264   fNtuFinalBestCandidates = new TNtuple("ntuFinalBestCandidates", "Final Best Candidates",  "run:event:muonTrack:nFinalCandidates:MCTrackRefExists:xVtx:yVtx:zVtx:motherPdg:triggerMatch:nClustersMC:nGoodClusters:pt:theta:eta:chi2AtPlane0:chi2AtPlane1:chi2AtPlane2:chi2AtPlane3:chi2AtPlane4:chi2AtPlane5:chi2AtPlane6:chi2AtPlane7:chi2AtPlane8:nClustersAtPlane0:nClustersAtPlane1:nClustersAtPlane2:nClustersAtPlane3:nClustersAtPlane4:nClustersAtPlane5:nClustersAtPlane6:nClustersAtPlane7:nClustersAtPlane8");
1265
1266 }
1267
1268 //============================================================================================================================================
1269
1270 void AliMuonForwardTrackFinder::SetTitleHistos() {
1271
1272   fHistRadiusEndOfAbsorber         -> SetXTitle("R_{abs}  [cm]");
1273   fHistNGoodClustersForFinalTracks -> SetXTitle("N_{GoodClusters}");
1274
1275   fHistDistanceGoodClusterFromTrackAtLastPlane                                  -> SetXTitle("Distance  [cm]");
1276   fHistDistanceGoodClusterFromTrackMinusDistanceBestClusterFromTrackAtLastPlane -> SetXTitle("Distance  [cm]");
1277
1278
1279   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
1280
1281     fHistNTracksAfterExtrapolation[iPlane] -> SetXTitle("N_{tracks}");
1282     fHistResearchRadius[iPlane]            -> SetXTitle("Research Radius  [cm]");
1283
1284     fHistChi2Cluster_GoodCluster[iPlane]         -> SetXTitle("#chi^{2}/ndf");
1285     fHistChi2Cluster_BadCluster[iPlane]          -> SetXTitle("#chi^{2}/ndf");
1286
1287     fHistGlobalChi2AtPlaneFor_GOOD_CandidatesOfTrackableMuons[iPlane] -> SetXTitle("#chi^{2}/ndf");
1288     fHistGlobalChi2AtPlaneFor_BAD_CandidatesOfTrackableMuons[iPlane]  -> SetXTitle("#chi^{2}/ndf");
1289
1290   }
1291
1292 }
1293
1294 //===========================================================================================================================================
1295
1296 void AliMuonForwardTrackFinder::BookPlanes() {
1297
1298   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
1299     fGrMFTPlane[kAllClusters][iPlane] = new TGraph();
1300     fGrMFTPlane[kAllClusters][iPlane] -> SetName(Form("fGrMFTPlane_%02d_AllClusters",iPlane));
1301     fGrMFTPlane[kAllClusters][iPlane] -> SetMarkerStyle(20);
1302     //    fGrMFTPlane[kAllClusters][iPlane] -> SetMarkerSize(0.5);
1303     //    fGrMFTPlane[kAllClusters][iPlane] -> SetMarkerSize(0.3);
1304     fGrMFTPlane[kAllClusters][iPlane] -> SetMarkerSize(0.2);
1305   }
1306
1307   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
1308     fGrMFTPlane[kClustersGoodChi2][iPlane] = new TGraph();
1309     fGrMFTPlane[kClustersGoodChi2][iPlane] -> SetName(Form("fGrMFTPlane_%02d_ClustersGoodChi2",iPlane));
1310     fGrMFTPlane[kClustersGoodChi2][iPlane] -> SetMarkerStyle(20);
1311     //    fGrMFTPlane[kClustersGoodChi2][iPlane] -> SetMarkerSize(0.8);
1312     //    fGrMFTPlane[kClustersGoodChi2][iPlane] -> SetMarkerSize(0.4);
1313     fGrMFTPlane[kClustersGoodChi2][iPlane] -> SetMarkerSize(0.3);
1314     fGrMFTPlane[kClustersGoodChi2][iPlane] -> SetMarkerColor(kBlue);
1315   }
1316
1317   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
1318     fGrMFTPlane[kClusterOfTrack][iPlane] = new TGraph();
1319     fGrMFTPlane[kClusterOfTrack][iPlane] -> SetName(Form("fGrMFTPlane_%02d_ClustersOfTrack",iPlane));
1320     fGrMFTPlane[kClusterOfTrack][iPlane] -> SetMarkerStyle(25);
1321     //    fGrMFTPlane[kClusterOfTrack][iPlane] -> SetMarkerSize(1.2);
1322     fGrMFTPlane[kClusterOfTrack][iPlane] -> SetMarkerSize(0.9);
1323     fGrMFTPlane[kClusterOfTrack][iPlane] -> SetMarkerColor(kRed);
1324     fGrMFTPlane[kClusterOfTrack][iPlane] -> SetTitle(Form("Plane %d (%3.1f cm)", iPlane, fZPlane[iPlane]));
1325   }
1326
1327   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
1328     fGrMFTPlane[kClusterCorrectMC][iPlane] = new TGraph();
1329     fGrMFTPlane[kClusterCorrectMC][iPlane] -> SetName(Form("fGrMFTPlane_%02d_ClustersCorrectMC",iPlane));
1330     fGrMFTPlane[kClusterCorrectMC][iPlane] -> SetMarkerStyle(20);
1331     //    fGrMFTPlane[kClusterCorrectMC][iPlane] -> SetMarkerSize(0.8);
1332     fGrMFTPlane[kClusterCorrectMC][iPlane] -> SetMarkerSize(0.5);
1333     fGrMFTPlane[kClusterCorrectMC][iPlane] -> SetMarkerColor(kGreen);
1334   }
1335
1336   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
1337     fCircleExt[iPlane] = new TEllipse(0., 0., fRPlaneMax[iPlane], fRPlaneMax[iPlane]);
1338     fCircleInt[iPlane] = new TEllipse(0., 0., fRPlaneMin[iPlane], fRPlaneMin[iPlane]);
1339   }
1340
1341   fTxtDummy = new TLatex(0.10, 0.59, "Best Candidate:");
1342
1343   //---------------------------------------------------
1344
1345   fMrkAllClust = new TMarker(0.10, 0.32, 20);
1346   fMrkAllClust -> SetMarkerSize(0.5);
1347
1348   fMrkClustGoodChi2 = new TMarker(0.10, 0.26, 20);
1349   fMrkClustGoodChi2 -> SetMarkerSize(0.8);
1350   fMrkClustGoodChi2 -> SetMarkerColor(kBlue);
1351
1352   fMrkClustMC = new TMarker(0.10, 0.20, 20);
1353   fMrkClustMC -> SetMarkerSize(0.8);
1354   fMrkClustMC -> SetMarkerColor(kGreen);
1355
1356   fMrkClustOfTrack = new TMarker(0.10, 0.14, 25);
1357   fMrkClustOfTrack -> SetMarkerSize(1.2);
1358   fMrkClustOfTrack -> SetMarkerColor(kRed);
1359
1360   fTxtAllClust = new TLatex(0.15, 0.30, "All Clusters");
1361   fTxtAllClust -> SetTextSize(0.040);
1362
1363   fTxtClustGoodChi2 = new TLatex(0.15, 0.24, "Clusters involved in the research");
1364   fTxtClustGoodChi2 -> SetTextSize(0.040);
1365
1366   fTxtClustMC = new TLatex(0.15, 0.18, "MC good clusters");
1367   fTxtClustMC -> SetTextSize(0.040);
1368
1369   fTxtClustOfTrack = new TLatex(0.15, 0.12, "Clusters of the best candidate");
1370   fTxtClustOfTrack -> SetTextSize(0.040);
1371
1372 }
1373
1374 //===========================================================================================================================================
1375
1376 void AliMuonForwardTrackFinder::ResetPlanes() {
1377
1378   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
1379     for (Int_t iGr=0; iGr<4; iGr++) {
1380       Int_t nOldClusters = fGrMFTPlane[iGr][iPlane]->GetN();
1381       for (Int_t iPoint=nOldClusters-1; iPoint>=0; iPoint--) fGrMFTPlane[iGr][iPlane]->RemovePoint(iPoint);
1382     }
1383   }
1384
1385 }
1386
1387 //===========================================================================================================================================
1388
1389 void AliMuonForwardTrackFinder::PrintParticleHistory() {
1390
1391   AliDebug(1, "Entering");
1392
1393   TString history = "";
1394
1395   TParticle *part = 0;
1396   if (fLabelMC>=0) part = fStack->Particle(fLabelMC);
1397
1398   AliDebug(1, Form("fStack->Particle(fLabelMC) = %p", part));
1399
1400   if (part) {
1401     if (part->GetFirstMother() != -1) {
1402       TParticle *partMother = fStack->Particle(part->GetFirstMother());
1403       AliDebug(1, Form("fStack->Particle(part->GetFirstMother() = %p", partMother));
1404       if (partMother) {
1405         Char_t newName[100];
1406         if (partMother->GetFirstMother() != -1) history += "...  #rightarrow ";
1407         PDGNameConverter(partMother->GetName(), newName);
1408         history += Form("%s #rightarrow ", newName);
1409       }
1410     }
1411     Char_t newName[100];
1412     PDGNameConverter(part->GetName(), newName);
1413     history += Form("%s  at  z = %5.1f cm", newName, part->Vz());
1414     //  printf("%s", history.Data());
1415   }
1416   else history += "NO AVAILABLE HISTORY";
1417
1418   fTxtMuonHistory = new TLatex(0.10, 0.86, history.Data());
1419
1420   // Filling particle history in the fFinalBestCandidate
1421
1422   if (part) {
1423     for (Int_t iParent=0; iParent<AliMuonForwardTrack::fgkNParentsMax; iParent++) {
1424       if (part->GetFirstMother() == -1) break;
1425       if (!(fStack->Particle(part->GetFirstMother()))) break;
1426       AliDebug(1, Form("fStack->Particle(part->GetFirstMother() = %p", fStack->Particle(part->GetFirstMother())));
1427       fFinalBestCandidate->SetParentMCLabel(iParent, part->GetFirstMother());
1428       fFinalBestCandidate->SetParentPDGCode(iParent, fStack->Particle(part->GetFirstMother())->GetPdgCode());
1429       part = fStack->Particle(part->GetFirstMother());
1430     }
1431   }
1432
1433 }
1434
1435 //===========================================================================================================================================
1436
1437 Bool_t AliMuonForwardTrackFinder::IsMother(const Char_t *nameMother) {
1438
1439   Bool_t result = kFALSE;
1440
1441   TParticle *part = 0;
1442   if (fLabelMC>=0) part = fStack->Particle(fLabelMC);
1443
1444   if (part) {
1445     if (part->GetFirstMother() != -1) {
1446       TParticle *partMother = fStack->Particle(part->GetFirstMother());
1447       if (partMother) {
1448         if (!strcmp(partMother->GetName(), nameMother)) result=kTRUE;
1449       }
1450     }
1451   }
1452
1453   return result;
1454
1455 }
1456
1457 //===========================================================================================================================================
1458
1459 void AliMuonForwardTrackFinder::DrawPlanes() {
1460
1461   fCanvas -> Clear();
1462   if (fNPlanesMFT <= 5)       fCanvas -> Divide(3,2);
1463   else if (fNPlanesMFT <= 11) fCanvas -> Divide(4,3);
1464   else if (fNPlanesMFT <= 19) fCanvas -> Divide(5,4);
1465
1466   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
1467
1468     fCanvas->cd(fNPlanesMFT-iPlane+1);
1469
1470     fGrMFTPlane[kClusterOfTrack][iPlane] -> GetXaxis() -> SetLimits(-1.1*fRPlaneMax[fNPlanesMFT-1], +1.1*fRPlaneMax[fNPlanesMFT-1]);
1471     fGrMFTPlane[kClusterOfTrack][iPlane] -> GetYaxis() -> SetRangeUser(-1.1*fRPlaneMax[fNPlanesMFT-1], +1.1*fRPlaneMax[fNPlanesMFT-1]);
1472     fGrMFTPlane[kClusterOfTrack][iPlane] -> GetXaxis() -> SetTitle("X  [cm]");
1473     fGrMFTPlane[kClusterOfTrack][iPlane] -> GetYaxis() -> SetTitle("Y  [cm]");
1474     fGrMFTPlane[kClusterOfTrack][iPlane] -> Draw("ap");
1475
1476     fCircleExt[iPlane] -> Draw("same");
1477     fCircleInt[iPlane] -> Draw("same");
1478
1479     if (fGrMFTPlane[kAllClusters][iPlane]->GetN())       fGrMFTPlane[kAllClusters][iPlane]      -> Draw("psame");
1480     if (fGrMFTPlane[kClustersGoodChi2][iPlane]->GetN())  fGrMFTPlane[kClustersGoodChi2][iPlane] -> Draw("psame");
1481     if (fGrMFTPlane[kClusterOfTrack][iPlane]->GetN())    fGrMFTPlane[kClusterOfTrack][iPlane]   -> Draw("psame");
1482     if (fGrMFTPlane[kClusterCorrectMC][iPlane]->GetN())  fGrMFTPlane[kClusterCorrectMC][iPlane] -> Draw("psame");
1483
1484     fTxtTrackChi2[iPlane] -> Draw("same");
1485
1486   }
1487
1488   fCanvas -> cd(1);
1489   fTxtMuonHistory       -> Draw();
1490   fTxtDummy             -> Draw("same");
1491   if (fMatchingMode==kRealMatching) fTxtTrackGoodClusters -> Draw("same");
1492   fTxtTrackFinalChi2    -> Draw("same");
1493   fTxtTrackMomentum     -> Draw("same");
1494   if (fMatchingMode==kRealMatching) fTxtFinalCandidates   -> Draw("same");
1495
1496   fMrkAllClust      -> Draw("same");
1497   fMrkClustGoodChi2 -> Draw("same");
1498   fMrkClustMC       -> Draw("same");
1499   fMrkClustOfTrack  -> Draw("same");
1500
1501   fTxtAllClust      -> Draw("same");
1502   fTxtClustGoodChi2 -> Draw("same");
1503   fTxtClustMC       -> Draw("same");
1504   fTxtClustOfTrack  -> Draw("same");
1505
1506   //  fCanvas -> SaveAs(Form("%s/figures/eventDisplay/run%d_event%d_track%d.eps", fOutDir.Data(), fRun, fEv, fCountRealTracksAnalyzedOfEvent));
1507   fCanvas -> SaveAs(Form("%s/figures/eventDisplay/run%d_event%d_track%d.gif", fOutDir.Data(), fRun, fEv, fCountRealTracksAnalyzedOfEvent));
1508   if (IsMother("phi")) {
1509     fCanvas -> SaveAs(Form("%s/figures/eventDisplay/run%d_event%d_track%d.phi.gif", fOutDir.Data(), fRun, fEv, fCountRealTracksAnalyzedOfEvent));
1510     fCanvas -> SaveAs(Form("%s/figures/eventDisplay/run%d_event%d_track%d.phi.eps", fOutDir.Data(), fRun, fEv, fCountRealTracksAnalyzedOfEvent));
1511   }
1512   if (IsMother("J/psi")) {
1513     fCanvas -> SaveAs(Form("%s/figures/eventDisplay/run%d_event%d_track%d.jPsi.gif", fOutDir.Data(), fRun, fEv, fCountRealTracksAnalyzedOfEvent));
1514     fCanvas -> SaveAs(Form("%s/figures/eventDisplay/run%d_event%d_track%d.jPsi.eps", fOutDir.Data(), fRun, fEv, fCountRealTracksAnalyzedOfEvent));
1515   }
1516
1517 }
1518
1519 //===========================================================================================================================================
1520
1521 void AliMuonForwardTrackFinder::Terminate() {
1522
1523   AliInfo("");
1524   AliInfo("---------------------------------------------------------------------------------------------------------------");
1525   AliInfo(Form("%8d  tracks analyzed",                                                     fCountRealTracksAnalyzed));
1526   AliInfo(Form("%8d  tracks with MC ref",                                                  fCountRealTracksWithRefMC));
1527   AliInfo(Form("%8d  tracks with MC ref & trigger match",                                  fCountRealTracksWithRefMC_andTrigger));
1528   if (fMatchingMode==kRealMatching) {
1529     AliInfo(Form("%8d  tracks analyzed with final candidates",                             fCountRealTracksAnalyzedWithFinalCandidates));
1530   }
1531   else {
1532     AliInfo(Form("%8d  tracks matched with their MC clusters",                             fCountRealTracksAnalyzedWithFinalCandidates));
1533   }
1534 //   printf("%8d  tracks with MC ref & trigger match & pt>%3.1f GeV/c",                 fCountRealTracksWithRefMC_andTrigger_andGoodPt, fLowPtCut);
1535 //   printf("%8d  tracks with MC ref & trigger match & pt>%3.1f GeV/c & correct R_abs", fCountRealTracksWithRefMC_andTrigger_andGoodPt_andGoodTheta, fLowPtCut);
1536   AliInfo("---------------------------------------------------------------------------------------------------------------");
1537
1538   WriteOutputTree();
1539   WriteHistos();
1540
1541 }
1542
1543 //==========================================================================================================================================
1544
1545 void AliMuonForwardTrackFinder::FillOutputTree() {
1546
1547   if (!fMuonForwardTracks || !fOutputEventTree) return;
1548
1549   AliDebug(1, Form("Filling output tree %p with %p having %d entries whose 1st entry is %p",
1550                    fOutputEventTree, fMuonForwardTracks, fMuonForwardTracks->GetEntries(), fMuonForwardTracks->At(0)));
1551
1552   //  fOutputTreeFile->cd();
1553   fOutputEventTree->Fill();
1554   AliDebug(1, Form("\nFilled Tree: nEvents = %d!!!!\n", Int_t(fOutputEventTree->GetEntries())));
1555
1556 }
1557
1558 //==========================================================================================================================================
1559
1560 void AliMuonForwardTrackFinder::WriteOutputTree() {
1561
1562   if (!fOutputEventTree || !fOutputTreeFile) return;
1563
1564   fOutputTreeFile -> cd();
1565
1566   fOutputEventTree -> Write();
1567   fOutputTreeFile -> Close();
1568
1569 }
1570
1571 //==========================================================================================================================================
1572
1573 void AliMuonForwardTrackFinder::WriteHistos() {
1574
1575   fOutputQAFile = new TFile(Form("MuonGlobalTracking.QA.run%d.root", fRun), "recreate");
1576   fOutputQAFile -> cd();
1577
1578   fHistRadiusEndOfAbsorber         -> Write();
1579   fHistNGoodClustersForFinalTracks -> Write();
1580
1581   fHistDistanceGoodClusterFromTrackAtLastPlane                                  -> Write();
1582   fHistDistanceGoodClusterFromTrackMinusDistanceBestClusterFromTrackAtLastPlane -> Write();
1583
1584   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
1585
1586     fHistNTracksAfterExtrapolation[iPlane] -> Write();
1587     fHistResearchRadius[iPlane]            -> Write();
1588
1589     fHistChi2Cluster_GoodCluster[iPlane]   -> Write();
1590     fHistChi2Cluster_BadCluster[iPlane]    -> Write();
1591
1592     fHistGlobalChi2AtPlaneFor_GOOD_CandidatesOfTrackableMuons[iPlane] -> Write();
1593     fHistGlobalChi2AtPlaneFor_BAD_CandidatesOfTrackableMuons[iPlane]  -> Write();
1594
1595   }
1596
1597   fNtuFinalCandidates     -> Write();
1598   fNtuFinalBestCandidates -> Write();
1599
1600   fOutputQAFile -> Close();
1601
1602 }
1603
1604 //===========================================================================================================================================
1605
1606 void AliMuonForwardTrackFinder::PDGNameConverter(const Char_t *nameIn, Char_t *nameOut) {
1607
1608   if      (!strcmp(nameIn, "mu+"))     snprintf(nameOut, 50, "#mu^{+}");
1609   else if (!strcmp(nameIn, "mu-"))     snprintf(nameOut, 50, "#mu^{-}");
1610   else if (!strcmp(nameIn, "pi+"))     snprintf(nameOut, 50, "#pi^{+}");
1611   else if (!strcmp(nameIn, "pi-"))     snprintf(nameOut, 50, "#pi^{-}");
1612   else if (!strcmp(nameIn, "K+"))      snprintf(nameOut, 50, "K^{+}");
1613   else if (!strcmp(nameIn, "K-"))      snprintf(nameOut, 50, "K^{-}");
1614   else if (!strcmp(nameIn, "K*+"))     snprintf(nameOut, 50, "K^{*+}");
1615   else if (!strcmp(nameIn, "K*-"))     snprintf(nameOut, 50, "K^{*-}");
1616   else if (!strcmp(nameIn, "K_S0"))    snprintf(nameOut, 50, "K_{S}^{0}");
1617   else if (!strcmp(nameIn, "K_L0"))    snprintf(nameOut, 50, "K_{L}^{0}");
1618   else if (!strcmp(nameIn, "K0"))      snprintf(nameOut, 50, "K^{0}");
1619   else if (!strcmp(nameIn, "K0_bar"))  snprintf(nameOut, 50, "#bar{K}^{0}");
1620   else if (!strcmp(nameIn, "K*0"))     snprintf(nameOut, 50, "K^{*0}");
1621   else if (!strcmp(nameIn, "K*0_bar")) snprintf(nameOut, 50, "#bar{K}^{*0}");
1622   else if (!strcmp(nameIn, "rho0"))    snprintf(nameOut, 50, "#rho^{0}");
1623   else if (!strcmp(nameIn, "rho+"))    snprintf(nameOut, 50, "#rho^{+}");
1624   else if (!strcmp(nameIn, "rho-"))    snprintf(nameOut, 50, "#rho^{-}");
1625   else if (!strcmp(nameIn, "omega"))   snprintf(nameOut, 50, "#omega");
1626   else if (!strcmp(nameIn, "eta'"))    snprintf(nameOut, 50, "#eta'");
1627   else if (!strcmp(nameIn, "phi"))     snprintf(nameOut, 50, "#phi");
1628
1629   else if (!strcmp(nameIn, "D-"))     snprintf(nameOut, 50, "D^{-}");
1630   else if (!strcmp(nameIn, "D+"))     snprintf(nameOut, 50, "D^{+}");
1631   else if (!strcmp(nameIn, "D0"))     snprintf(nameOut, 50, "D^{0}");
1632   else if (!strcmp(nameIn, "D0_bar")) snprintf(nameOut, 50, "#bar{D}^{0}");
1633   else if (!strcmp(nameIn, "D*-"))    snprintf(nameOut, 50, "D^{*-}");
1634   else if (!strcmp(nameIn, "D*+"))    snprintf(nameOut, 50, "D^{*+}");
1635   else if (!strcmp(nameIn, "D_s+"))   snprintf(nameOut, 50, "D_{s}^{+}");
1636   else if (!strcmp(nameIn, "D*_s+"))  snprintf(nameOut, 50, "D_{s}^{*+}");
1637
1638   else if (!strcmp(nameIn, "B-"))       snprintf(nameOut, 50, "B^{-}");
1639   else if (!strcmp(nameIn, "B+"))       snprintf(nameOut, 50, "B^{+}");
1640   else if (!strcmp(nameIn, "B_s0_bar")) snprintf(nameOut, 50, "#bar{B}_{s}^{0}");
1641
1642   else if (!strcmp(nameIn, "antiproton"))  snprintf(nameOut, 50, "#bar{p}");
1643   else if (!strcmp(nameIn, "proton"))      snprintf(nameOut, 50, "p");
1644   else if (!strcmp(nameIn, "neutron"))     snprintf(nameOut, 50, "n");
1645   else if (!strcmp(nameIn, "Sigma+"))      snprintf(nameOut, 50, "#Sigma^{+}");
1646   else if (!strcmp(nameIn, "Delta+"))      snprintf(nameOut, 50, "#Delta{+}");
1647   else if (!strcmp(nameIn, "Delta--"))     snprintf(nameOut, 50, "#Delta{--}");
1648   else if (!strcmp(nameIn, "Lambda0"))     snprintf(nameOut, 50, "#Lambda_0");
1649   else if (!strcmp(nameIn, "Lambda0_bar")) snprintf(nameOut, 50, "#bar{Lambda}_0");
1650
1651   else snprintf(nameOut, 50, "%s", nameIn);
1652
1653 }
1654
1655 //===========================================================================================================================================
1656
1657 void AliMuonForwardTrackFinder::SetDraw(Bool_t drawOption) {
1658
1659   fDrawOption = drawOption;
1660
1661   if (!fCanvas) {
1662     fCanvas = new TCanvas("tracking", "tracking", 1200, 800);
1663     fCanvas -> Divide(3,2);
1664   }
1665
1666 }
1667
1668 //===========================================================================================================================================
1669
1670 Bool_t AliMuonForwardTrackFinder::InitGRP() {
1671
1672   //------------------------------------
1673   // Initialization of the GRP entry
1674   //------------------------------------
1675
1676   AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/GRP/Data");
1677
1678   if (entry) {
1679
1680     TMap* m = dynamic_cast<TMap*>(entry->GetObject());  // old GRP entry
1681
1682     if (m) {
1683        AliInfo("Found a TMap in GRP/GRP/Data, converting it into an AliGRPObject");
1684        m->Print();
1685        fGRPData = new AliGRPObject();
1686        fGRPData->ReadValuesFromMap(m);
1687     }
1688
1689     else {
1690        AliInfo("Found an AliGRPObject in GRP/GRP/Data, reading it");
1691        fGRPData = dynamic_cast<AliGRPObject*>(entry->GetObject());  // new GRP entry
1692        entry->SetOwner(0);
1693     }
1694
1695     //    FIX ME: The unloading of GRP entry is temporarily disabled
1696     //    because ZDC and VZERO are using it in order to initialize
1697     //    their reconstructor objects. In the future one has to think
1698     //    of propagating AliRunInfo to the reconstructors.
1699     //    AliCDBManager::Instance()->UnloadFromCache("GRP/GRP/Data");
1700   }
1701
1702   if (!fGRPData) {
1703      AliError("No GRP entry found in OCDB!");
1704      return kFALSE;
1705   }
1706
1707   TString lhcState = fGRPData->GetLHCState();
1708   if (lhcState==AliGRPObject::GetInvalidString()) {
1709     AliError("GRP/GRP/Data entry:  missing value for the LHC state ! Using UNKNOWN");
1710     lhcState = "UNKNOWN";
1711   }
1712
1713   TString beamType = fGRPData->GetBeamType();
1714   if (beamType==AliGRPObject::GetInvalidString()) {
1715     AliError("GRP/GRP/Data entry:  missing value for the beam type ! Using UNKNOWN");
1716     beamType = "UNKNOWN";
1717   }
1718
1719   Float_t beamEnergy = fGRPData->GetBeamEnergy();
1720   if (beamEnergy==AliGRPObject::GetInvalidFloat()) {
1721     AliError("GRP/GRP/Data entry:  missing value for the beam energy ! Using 0");
1722     beamEnergy = 0;
1723   }
1724
1725   TString runType = fGRPData->GetRunType();
1726   if (runType==AliGRPObject::GetInvalidString()) {
1727     AliError("GRP/GRP/Data entry:  missing value for the run type ! Using UNKNOWN");
1728     runType = "UNKNOWN";
1729   }
1730
1731   Int_t activeDetectors = fGRPData->GetDetectorMask();
1732   if (activeDetectors==AliGRPObject::GetInvalidUInt()) {
1733     AliError("GRP/GRP/Data entry:  missing value for the detector mask ! Using 1074790399");
1734     activeDetectors = 1074790399;
1735   }
1736   AliDebug(1, Form("activeDetectors = %d", activeDetectors));
1737
1738   fRunInfo = new AliRunInfo(lhcState, beamType, beamEnergy, runType, activeDetectors);
1739   fRunInfo->Dump();
1740
1741   // *** Dealing with the magnetic field map
1742
1743   if ( TGeoGlobalMagField::Instance()->IsLocked() ) {
1744     if (TGeoGlobalMagField::Instance()->GetField()->TestBit(AliMagF::kOverrideGRP)) {
1745       AliInfo("ExpertMode!!! GRP information will be ignored !");
1746       AliInfo("ExpertMode!!! Running with the externally locked B field !");
1747     }
1748     else {
1749       AliInfo("Destroying existing B field instance!");
1750       delete TGeoGlobalMagField::Instance();
1751     }
1752   }
1753   if ( !TGeoGlobalMagField::Instance()->IsLocked() ) {
1754     // Construct the field map out of the information retrieved from GRP.
1755     Bool_t ok = kTRUE;
1756     // L3
1757     Float_t l3Current = fGRPData->GetL3Current((AliGRPObject::Stats)0);
1758     if (l3Current == AliGRPObject::GetInvalidFloat()) {
1759       AliError("GRP/GRP/Data entry:  missing value for the L3 current !");
1760       ok = kFALSE;
1761     }
1762
1763     Char_t l3Polarity = fGRPData->GetL3Polarity();
1764     if (l3Polarity == AliGRPObject::GetInvalidChar()) {
1765       AliError("GRP/GRP/Data entry:  missing value for the L3 polarity !");
1766       ok = kFALSE;
1767     }
1768
1769     // Dipole
1770     Float_t diCurrent = fGRPData->GetDipoleCurrent((AliGRPObject::Stats)0);
1771     if (diCurrent == AliGRPObject::GetInvalidFloat()) {
1772       AliError("GRP/GRP/Data entry:  missing value for the dipole current !");
1773       ok = kFALSE;
1774     }
1775
1776     Char_t diPolarity = fGRPData->GetDipolePolarity();
1777     if (diPolarity == AliGRPObject::GetInvalidChar()) {
1778       AliError("GRP/GRP/Data entry:  missing value for the dipole polarity !");
1779       ok = kFALSE;
1780     }
1781
1782     // read special bits for the polarity convention and map type
1783     Int_t  polConvention = fGRPData->IsPolarityConventionLHC() ? AliMagF::kConvLHC : AliMagF::kConvDCS2008;
1784     Bool_t uniformB = fGRPData->IsUniformBMap();
1785
1786     if (ok) {
1787       AliMagF* fld = AliMagF::CreateFieldMap(TMath::Abs(l3Current) * (l3Polarity ? -1:1),
1788                                              TMath::Abs(diCurrent) * (diPolarity ? -1:1),
1789                                              polConvention,uniformB,beamEnergy, beamType.Data());
1790       if (fld) {
1791         TGeoGlobalMagField::Instance()->SetField( fld );
1792         TGeoGlobalMagField::Instance()->Lock();
1793         AliInfo("Running with the B field constructed out of GRP !");
1794       }
1795       else AliFatal("Failed to create a B field map !");
1796     }
1797     else AliFatal("B field is neither set nor constructed from GRP ! Exitig...");
1798   }
1799
1800   return kTRUE;
1801 }
1802
1803 //====================================================================================================================================================
1804
1805 Bool_t AliMuonForwardTrackFinder::SetRunNumber() {
1806
1807   AliCDBManager *man = AliCDBManager::Instance();
1808
1809   if (!fRunLoader) {
1810     AliError("No run loader found!");
1811     return kFALSE;
1812   }
1813   else {
1814     fRunLoader->LoadHeader();
1815     // read run number from gAlice
1816     if (fRunLoader->GetHeader()) {
1817       man->SetRun(fRunLoader->GetHeader()->GetRun());
1818       fRunLoader->UnloadHeader();
1819     }
1820     else {
1821       AliError("No run-loader header found!");
1822       return kFALSE;
1823     }
1824   }
1825
1826   return kTRUE;
1827
1828 }
1829
1830 //====================================================================================================================================================
1831