MFT/AliMFTTrackerMU.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 //                          MFT tracker
  19 //
  20 // Class for the creation of the "global muon tracks" built from the tracks reconstructed in the
  21 // muon spectrometer and the clusters of the Muon Forward Tracker
  22 //
  23 // Contact author: antonio.uras@cern.ch
  24 //
  25 //====================================================================================================================================================
  26
  27 #include "TTree.h"
  28 #include "AliLog.h"
  29 #include "AliGeomManager.h"
  30 #include "AliESDEvent.h"
  31 #include "AliESDMuonTrack.h"
  32 #include "AliESDMuonGlobalTrack.h"
  33 #include "AliMFTTrackerMU.h"
  34 #include "TMath.h"
  35 #include "AliRun.h"
  36 #include "AliRunLoader.h"
  37 #include "AliHeader.h"
  38 #include "AliGenEventHeader.h"
  39 #include "TArrayF.h"
  40 #include "AliMFT.h"
  41 #include "AliMUONTrackExtrap.h"
  42 #include "AliMUONTrack.h"
  43 #include "AliMUONESDInterface.h"
  44 #include "AliMuonForwardTrack.h"
  45 #include "AliMUONConstants.h"
  46
  47 ClassImp(AliMFTTrackerMU)
  48
  49 const Double_t AliMFTTrackerMU::fRadLengthSi = AliMFTConstants::fRadLengthSi;
  50
  51 //====================================================================================================================================================
  52
  53 AliMFTTrackerMU::AliMFTTrackerMU() :
  54   AliTracker(),
  55   fESD(0),
  56   fMFT(0),
  57   fSegmentation(0),
  58   fNPlanesMFT(0),
  59   fNPlanesMFTAnalyzed(0),
  60   fSigmaClusterCut(2),
  61   fScaleSigmaClusterCut(1.),
  62   fNMaxMissingMFTClusters(0),
  63   fGlobalTrackingDiverged(kFALSE),
  64   fCandidateTracks(0),
  65   fMUONTrack(0),
  66   fCurrentTrack(0),
  67   fFinalBestCandidate(0),
  68   fXExtrapVertex(0),
  69   fYExtrapVertex(0),
  70   fZExtrapVertex(0),
  71   fXExtrapVertexError(0),
  72   fYExtrapVertexError(0),
  73   fBransonCorrection(kFALSE)
  74 {
  75
  76   //--------------------------------------------------------------------
  77   // This is the AliMFTTrackerMU constructor
  78   //--------------------------------------------------------------------
  79
  80   fMFT = (AliMFT*) gAlice->GetDetector("MFT");
  81   fSegmentation = fMFT->GetSegmentation();
  82   SetNPlanesMFT(fSegmentation->GetNPlanes());
  83   AliMUONTrackExtrap::SetField();                 // set the magnetic field for track extrapolations
  84
  85   AliInfo(Form("fMFT = %p, fSegmentation = %p", fMFT, fSegmentation));
  86
  87   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
  88     fMFTClusterArray[iPlane]      = new TClonesArray("AliMFTCluster");
  89     fMFTClusterArrayFront[iPlane] = new TClonesArray("AliMFTCluster");
  90     fMFTClusterArrayBack[iPlane]  = new TClonesArray("AliMFTCluster");
  91     fMFTClusterArray[iPlane]      -> SetOwner(kTRUE);
  92     fMFTClusterArrayFront[iPlane] -> SetOwner(kTRUE);
  93     fMFTClusterArrayBack[iPlane]  -> SetOwner(kTRUE);
  94     fMinResearchRadiusAtPlane[iPlane] = 0.;
  95   }
  96
  97   fCandidateTracks = new TClonesArray("AliMuonForwardTrack",50000);
  98
  99 }
 100
 101 //====================================================================================================================================================
 102
 103 AliMFTTrackerMU::~AliMFTTrackerMU() {
 104
 105   // destructor
 106
 107   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
 108     fMFTClusterArray[iPlane] -> Delete();
 109     delete fMFTClusterArray[iPlane];
 110     delete fMFTClusterArrayFront[iPlane];
 111     delete fMFTClusterArrayBack[iPlane];
 112   }
 113
 114   delete fCandidateTracks;
 115
 116 }
 117
 118 //====================================================================================================================================================
 119
 120 Int_t AliMFTTrackerMU::LoadClusters(TTree *cTree) {
 121
 122   //--------------------------------------------------------------------
 123   // This function loads the MFT clusters
 124   //--------------------------------------------------------------------
 125
 126   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
 127     AliDebug(1, Form("Setting Address for Branch Plane_%02d", iPlane));
 128     cTree->SetBranchAddress(Form("Plane_%02d",iPlane), &fMFTClusterArray[iPlane]);
 129   }
 130
 131   if (!cTree->GetEvent()) return kFALSE;
 132   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
 133     AliInfo(Form("plane %02d: nClusters = %d\n", iPlane, fMFTClusterArray[iPlane]->GetEntries()));
 134   }
 135   SeparateFrontBackClusters();
 136
 137   return 0;
 138
 139 }
 140
 141 //====================================================================================================================================================
 142
 143 void AliMFTTrackerMU::UnloadClusters() {
 144
 145   //--------------------------------------------------------------------
 146   // This function unloads MFT clusters
 147   //--------------------------------------------------------------------
 148
 149   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
 150     fMFTClusterArray[iPlane]      -> Clear("C");
 151     fMFTClusterArrayFront[iPlane] -> Clear("C");
 152     fMFTClusterArrayBack[iPlane]  -> Clear("C");
 153   }
 154
 155 }
 156
 157 //====================================================================================================================================================
 158
 159 Int_t AliMFTTrackerMU::Clusters2Tracks(AliESDEvent *event) {
 160
 161   //--------------------------------------------------------------------
 162   // This functions reconstructs the Muon Forward Tracks
 163   // The clusters must be already loaded !
 164   //--------------------------------------------------------------------
 165
 166   // Tree of AliMuonForwardTrack objects. Created outside the ESD framework for cross-check purposes
 167
 168   TFile *outputFileMuonGlobalTracks = new TFile("MuonGlobalTracks.root", "update");
 169   TTree *outputTreeMuonGlobalTracks = new TTree("AliMuonForwardTracks", "Tree of AliMuonForwardTracks");
 170   TClonesArray *muonForwardTracks = new TClonesArray("AliMuonForwardTrack");
 171   outputTreeMuonGlobalTracks -> Branch("tracks", &muonForwardTracks);
 172
 173   //--------------------------------------------------------------------
 174
 175   fESD = event;
 176
 177   fXExtrapVertex = 0;
 178   fYExtrapVertex = 0;
 179   fZExtrapVertex = 0;
 180   fXExtrapVertexError = AliMFTConstants::fXVertexTolerance;
 181   fYExtrapVertexError = AliMFTConstants::fYVertexTolerance;
 182
 183   // Imposing a fixed primary vertex, in order to keep memory of its position. Taking the primary vertex would imply the risk
 184   // of loosing memory of its position when passing from ESD to AOD, due to possible refitting
 185   const AliESDVertex* esdVert = fESD->GetVertex();
 186   if (esdVert->GetNContributors() > 0 || !strcmp(esdVert->GetTitle(),"vertexer: smearMC")) {
 187     fXExtrapVertex = esdVert->GetX();
 188     fYExtrapVertex = esdVert->GetY();
 189     fZExtrapVertex = esdVert->GetZ();
 190     fXExtrapVertexError = TMath::Max(AliMFTConstants::fXVertexTolerance, esdVert->GetXRes());
 191     fYExtrapVertexError = TMath::Max(AliMFTConstants::fYVertexTolerance, esdVert->GetYRes());
 192     AliInfo(Form("Found ESD vertex from %d contributors (%f +/- %f,  %f +/- %f,  %f)",
 193                  esdVert->GetNContributors(),fXExtrapVertex,fXExtrapVertexError,fYExtrapVertex,fYExtrapVertexError,fZExtrapVertex));
 194   }
 195   else GetVertexFromMC();
 196
 197   //----------- Read ESD MUON tracks -------------------
 198
 199   Int_t nTracksMUON = event->GetNumberOfMuonTracks();
 200
 201   AliInfo(Form("Number of ESD MUON tracks: %d\n", nTracksMUON));
 202
 203   Int_t iTrack=0;
 204   while (iTrack<nTracksMUON) {
 205
 206     fNPlanesMFTAnalyzed = 0;
 207
 208     AliInfo("****************************************************************************************");
 209     AliInfo(Form("***************************   MUON TRACK %3d/%d   ***************************************", iTrack, nTracksMUON));
 210     AliInfo("****************************************************************************************");
 211
 212     fCandidateTracks -> Delete();
 213
 214     fNPlanesMFTAnalyzed = 0;
 215
 216     const AliESDMuonTrack *esdTrack = event->GetMuonTrack(iTrack);
 217     if (fMUONTrack) delete fMUONTrack;
 218     fMUONTrack = new AliMUONTrack();
 219
 220     AliMUONESDInterface::ESDToMUON(*esdTrack, *fMUONTrack, kFALSE);
 221
 222     if (!fMUONTrack->GetTrackParamAtCluster()->First()) {
 223       AliInfo("Skipping track, no parameters available!!!");
 224       iTrack++;
 225       continue;
 226     }
 227
 228     // the track we are going to build, starting from fMUONTrack and adding the MFT clusters
 229     AliMuonForwardTrack *track = new ((*fCandidateTracks)[0]) AliMuonForwardTrack();
 230     track -> SetMUONTrack(new AliMUONTrack(*fMUONTrack));
 231     track -> SetMCLabel(fMUONTrack->GetMCLabel());
 232     track -> SetMatchTrigger(fMUONTrack->GetMatchTrigger());
 233
 234     // track parameters linearly extrapolated from the first tracking station to the end of the absorber
 235     AliMUONTrackParam trackParamEndOfAbsorber(*((AliMUONTrackParam*)(fMUONTrack->GetTrackParamAtCluster()->First())));
 236     AliMUONTrackExtrap::ExtrapToZCov(&trackParamEndOfAbsorber, AliMUONConstants::AbsZEnd());   // absorber extends from -90 to -503 cm
 237     Double_t xEndOfAbsorber = trackParamEndOfAbsorber.GetNonBendingCoor();
 238     Double_t yEndOfAbsorber = trackParamEndOfAbsorber.GetBendingCoor();
 239     Double_t rAbsorber      = TMath::Sqrt(xEndOfAbsorber*xEndOfAbsorber + yEndOfAbsorber*yEndOfAbsorber);
 240     track -> SetRAtAbsorberEnd(rAbsorber);
 241
 242     //------------------------- NOW THE CYCLE OVER THE MFT PLANES STARTS ---------------------------------------
 243
 244     for (Int_t iPlane=fNPlanesMFT-1; iPlane>=0; iPlane--) {   /* *** do not reverse the order of this cycle!!!
 245                                                                  *** this reflects the fact that the extrapolation is performed
 246                                                                  *** starting from the last MFT plane back to the origin */
 247
 248       // --------- updating the array of candidates according to the clusters available in the i-th plane ---------
 249
 250       fNPlanesMFTAnalyzed++;
 251
 252       Int_t nCandidates = fCandidateTracks->GetEntriesFast();
 253       for (Int_t iCandidate=0; iCandidate<nCandidates; iCandidate++) {
 254         fCurrentTrack = (AliMuonForwardTrack*) fCandidateTracks->UncheckedAt(iCandidate);
 255
 256         // if the old track is compatible with the new cluster, the track is updated and inserted as new track in the array
 257         // (several new tracks can be created for one old track)
 258         if (FindClusterInPlane(iPlane) == kDiverged) {
 259           fGlobalTrackingDiverged = kTRUE;
 260           break;
 261         }
 262
 263         if ((fNPlanesMFTAnalyzed-fCurrentTrack->GetNMFTClusters())>fNMaxMissingMFTClusters || fIsPlaneMandatory[iPlane]) {
 264           fCandidateTracks->Remove(fCurrentTrack);     // the old track is removed after the check;
 265         }
 266       }
 267       if (fGlobalTrackingDiverged) {
 268         if (fScaleSigmaClusterCut>0) fScaleSigmaClusterCut -= 0.1;
 269         continue;
 270       }
 271
 272       fCandidateTracks->Compress();
 273
 274     }
 275
 276     // -------------------------- END OF THE CYCLE OVER THE MFT PLANES --------------------------------------------
 277
 278     fGlobalTrackingDiverged = kFALSE;
 279     fScaleSigmaClusterCut = 1.0;
 280
 281     AliDebug(1, "Finished cycle over planes");
 282
 283     iTrack++;
 284
 285     // If we have several final tracks, we must find the best candidate:
 286
 287     Int_t nFinalTracks = fCandidateTracks->GetEntriesFast();
 288     AliInfo(Form("nFinalTracks = %d", nFinalTracks));
 289
 290     Int_t nGoodClustersBestCandidate =  0;
 291     Int_t idBestCandidate            =  0;
 292     Double_t bestChi2                = -1.;  // variable defining the best candidate
 293
 294     for (Int_t iFinalCandidate=0; iFinalCandidate<nFinalTracks; iFinalCandidate++) {
 295
 296       AliMuonForwardTrack *finalTrack = (AliMuonForwardTrack*) fCandidateTracks->UncheckedAt(iFinalCandidate);
 297       Int_t nMFTClusters  = finalTrack->GetNMFTClusters();
 298
 299       Double_t chi2 = 0;
 300       for (Int_t iCluster=0; iCluster<nMFTClusters; iCluster++) {
 301         AliMFTCluster *localCluster = finalTrack->GetMFTCluster(iCluster);
 302         chi2 += localCluster->GetLocalChi2();
 303       }
 304       chi2 /= nMFTClusters;
 305
 306       // now comparing the tracks in order to find the best one
 307
 308       if (chi2<bestChi2 || bestChi2<0) {
 309         bestChi2 = chi2;
 310         idBestCandidate = iFinalCandidate;
 311       }
 312
 313     }
 314
 315     if (nFinalTracks) {
 316
 317       AliMuonForwardTrack *newTrack = (AliMuonForwardTrack*) fCandidateTracks->UncheckedAt(idBestCandidate);
 318       newTrack -> SetNWrongClustersMC(newTrack->GetNMFTClusters() - nGoodClustersBestCandidate);
 319
 320       new ((*muonForwardTracks)[muonForwardTracks->GetEntries()]) AliMuonForwardTrack(*newTrack);
 321
 322       //----------------------- Save the information to the AliESDMuonGlobalTrack object
 323
 324       newTrack -> EvalKinem(AliMFTConstants::fZEvalKinem);
 325
 326       AliDebug(2,"Creating a new Muon Global Track");
 327       AliESDMuonGlobalTrack *myESDTrack = event->NewMuonGlobalTrack();
 328       myESDTrack -> SetPxPyPz(newTrack->Px(), newTrack->Py(), newTrack->Pz());
 329
 330       myESDTrack -> SetLabel(newTrack->GetMCLabel());
 331       myESDTrack -> SetChi2OverNdf(newTrack->GetChi2OverNdf());
 332       myESDTrack -> SetCharge(newTrack->GetCharge());
 333       myESDTrack -> SetMatchTrigger(newTrack->GetMatchTrigger());
 334       myESDTrack -> SetNMFTClusters(newTrack->GetNMFTClusters());
 335       myESDTrack -> SetNWrongMFTClustersMC(newTrack->GetNWrongClustersMC());
 336       myESDTrack -> SetFirstTrackingPoint(newTrack->GetMFTCluster(0)->GetX(), newTrack->GetMFTCluster(0)->GetY(), newTrack->GetMFTCluster(0)->GetZ());
 337       myESDTrack -> SetXYAtVertex(newTrack->GetOffsetX(0., AliMFTConstants::fZEvalKinem), newTrack->GetOffsetY(0., AliMFTConstants::fZEvalKinem));
 338       myESDTrack -> SetRAtAbsorberEnd(newTrack->GetRAtAbsorberEnd());
 339       myESDTrack -> SetCovariances(newTrack->GetTrackParamAtMFTCluster(0)->GetCovariances());
 340       myESDTrack -> SetChi2MatchTrigger(esdTrack->GetChi2MatchTrigger());
 341       myESDTrack -> SetMuonClusterMap(esdTrack->GetMuonClusterMap());
 342       myESDTrack -> SetHitsPatternInTrigCh(esdTrack->GetHitsPatternInTrigCh());
 343       myESDTrack -> SetHitsPatternInTrigChTrk(esdTrack->GetHitsPatternInTrigChTrk());
 344       myESDTrack -> Connected(esdTrack->IsConnected());
 345
 346       ULong_t mftClusterPattern = 0;
 347       for (Int_t iCluster=0; iCluster<newTrack->GetNMFTClusters(); iCluster++) {
 348         AliMFTCluster *localCluster = newTrack->GetMFTCluster(iCluster);
 349         mftClusterPattern |= 1 << localCluster->GetPlane();
 350         mftClusterPattern |= IsCorrectMatch(localCluster, newTrack->GetMCLabel()) << (fNMaxPlanes + localCluster->GetPlane());
 351       }
 352       myESDTrack -> SetMFTClusterPattern(mftClusterPattern);
 353
 354       //---------------------------------------------------------------------------------
 355
 356     }
 357
 358     fFinalBestCandidate = NULL;
 359
 360   }
 361
 362   outputTreeMuonGlobalTracks -> Fill();
 363
 364   Int_t myEventID = 0;
 365   while (outputFileMuonGlobalTracks->cd(Form("Event%d",myEventID))) myEventID++;
 366   outputFileMuonGlobalTracks -> mkdir(Form("Event%d",myEventID));
 367   outputFileMuonGlobalTracks -> cd(Form("Event%d",myEventID));
 368   outputTreeMuonGlobalTracks -> Write();
 369   outputFileMuonGlobalTracks -> Close();
 370
 371   muonForwardTracks -> Delete();
 372   delete muonForwardTracks;
 373
 374   return 0;
 375
 376 }
 377
 378 //=========================================================================================================================================
 379
 380 void AliMFTTrackerMU::SeparateFrontBackClusters() {
 381
 382   for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
 383     fMFTClusterArrayFront[iPlane]->Delete();
 384     fMFTClusterArrayBack[iPlane] ->Delete();
 385     for (Int_t iCluster=0; iCluster<fMFTClusterArray[iPlane]->GetEntries(); iCluster++) {
 386       AliMFTCluster *cluster = (AliMFTCluster*) fMFTClusterArray[iPlane]->At(iCluster);
 387       if (TMath::Abs(cluster->GetZ())<TMath::Abs(fSegmentation->GetPlane(iPlane)->GetZCenter())) {
 388         new ((*fMFTClusterArrayFront[iPlane])[fMFTClusterArrayFront[iPlane]->GetEntries()]) AliMFTCluster(*cluster);
 389       }
 390       else {
 391         new ((*fMFTClusterArrayBack[iPlane])[fMFTClusterArrayBack[iPlane]->GetEntries()]) AliMFTCluster(*cluster);
 392       }
 393     }
 394   }
 395
 396 }
 397
 398 //==========================================================================================================================================
 399
 400 Int_t AliMFTTrackerMU::FindClusterInPlane(Int_t planeId) {
 401
 402   // !!!!!!!!! coordinates and errors on the interaction vertex should be taken from the event itself (ITS) if available
 403
 404   // propagate track to plane #planeId (both to front and back active sensors)
 405   // look for compatible clusters
 406   // update TrackParam at found cluster (if any) using Kalman Filter
 407
 408   AliMUONTrackParam currentParamFront, currentParamBack, currentParamForResearchFront, currentParamForResearchBack;
 409
 410   if (planeId == fNPlanesMFT-1) {      // last plane of the telecope
 411     currentParamFront = (*((AliMUONTrackParam*)(fMUONTrack->GetTrackParamAtCluster()->First())));
 412     currentParamBack  = (*((AliMUONTrackParam*)(fMUONTrack->GetTrackParamAtCluster()->First())));
 413     currentParamForResearchFront = currentParamFront;
 414     currentParamForResearchBack  = currentParamBack;
 415     if (fBransonCorrection) {
 416       AliMUONTrackExtrap::ExtrapToVertex(&currentParamFront, fXExtrapVertex, fYExtrapVertex, fZExtrapVertex, fXExtrapVertexError, fYExtrapVertexError);
 417       AliMUONTrackExtrap::ExtrapToVertex(&currentParamBack,  fXExtrapVertex, fYExtrapVertex, fZExtrapVertex, fXExtrapVertexError, fYExtrapVertexError);
 418     }
 419     else {
 420       AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(&currentParamFront, fZExtrapVertex);
 421       AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(&currentParamBack,  fZExtrapVertex);
 422     }
 423     AliMUONTrackExtrap::ExtrapToVertex(&currentParamForResearchFront, fXExtrapVertex, fYExtrapVertex, fZExtrapVertex, fXExtrapVertexError, fYExtrapVertexError);
 424     AliMUONTrackExtrap::ExtrapToVertex(&currentParamForResearchBack,  fXExtrapVertex, fYExtrapVertex, fZExtrapVertex, fXExtrapVertexError, fYExtrapVertexError);
 425   }
 426   else {          // MFT planes others than the last one: mult. scattering correction because of the upstream MFT planes is performed
 427     currentParamFront = (*((AliMUONTrackParam*)(fCurrentTrack->GetTrackParamAtCluster()->First())));
 428     currentParamBack  = (*((AliMUONTrackParam*)(fCurrentTrack->GetTrackParamAtCluster()->First())));
 429     currentParamForResearchFront = currentParamFront;
 430     currentParamForResearchBack  = currentParamBack;
 431     AliMUONTrackExtrap::AddMCSEffect(&currentParamFront,           (fSegmentation->GetPlane(planeId+1)->GetEquivalentSilicon()+
 432                                                                     fSegmentation->GetPlane(planeId)->GetEquivalentSiliconBeforeFront())/fRadLengthSi,-1.);
 433     AliMUONTrackExtrap::AddMCSEffect(&currentParamForResearchFront,(fSegmentation->GetPlane(planeId+1)->GetEquivalentSilicon()+
 434                                                                     fSegmentation->GetPlane(planeId)->GetEquivalentSiliconBeforeFront())/fRadLengthSi,-1.);
 435     AliMUONTrackExtrap::AddMCSEffect(&currentParamBack,            (fSegmentation->GetPlane(planeId+1)->GetEquivalentSilicon()+
 436                                                                     fSegmentation->GetPlane(planeId)->GetEquivalentSiliconBeforeBack())/fRadLengthSi,-1.);
 437     AliMUONTrackExtrap::AddMCSEffect(&currentParamForResearchBack, (fSegmentation->GetPlane(planeId+1)->GetEquivalentSilicon()+
 438                                                                     fSegmentation->GetPlane(planeId)->GetEquivalentSiliconBeforeBack())/fRadLengthSi,-1.);
 439   }
 440   // for all planes: extrapolation to the Z of the plane
 441   AliMUONTrackExtrap::ExtrapToZCov(&currentParamFront,            -1.*fSegmentation->GetPlane(planeId)->GetZCenterActiveFront());
 442   AliMUONTrackExtrap::ExtrapToZCov(&currentParamForResearchFront, -1.*fSegmentation->GetPlane(planeId)->GetZCenterActiveFront());
 443   AliMUONTrackExtrap::ExtrapToZCov(&currentParamBack,             -1.*fSegmentation->GetPlane(planeId)->GetZCenterActiveBack());
 444   AliMUONTrackExtrap::ExtrapToZCov(&currentParamForResearchBack,  -1.*fSegmentation->GetPlane(planeId)->GetZCenterActiveBack());
 445
 446   //---------------------------------------------------------------------------------------
 447
 448   TMatrixD covFront(5,5); covFront = currentParamForResearchFront.GetCovariances();
 449   TMatrixD covBack(5,5);  covBack  = currentParamForResearchBack.GetCovariances();
 450
 451   Double_t squaredError_X_Front = covFront(0,0);
 452   Double_t squaredError_Y_Front = covFront(2,2);
 453   Double_t squaredError_X_Back  = covBack(0,0);
 454   Double_t squaredError_Y_Back  = covBack(2,2);
 455
 456   Double_t corrFact = 1.0;
 457
 458   Double_t researchRadiusFront = TMath::Sqrt(squaredError_X_Front + squaredError_Y_Front);
 459   Double_t researchRadiusBack  = TMath::Sqrt(squaredError_X_Back  + squaredError_Y_Back);
 460   if (0.5*(researchRadiusFront+researchRadiusBack)<fMinResearchRadiusAtPlane[planeId]) {
 461     corrFact = fMinResearchRadiusAtPlane[planeId]/(0.5*(researchRadiusFront+researchRadiusBack));
 462   }
 463
 464   //---------------------------------------------------------------------------------------
 465
 466   Double_t chi2cut = 2.*fScaleSigmaClusterCut*fScaleSigmaClusterCut*fSigmaClusterCut*fSigmaClusterCut;     // depends on the number of variables (here, 2)
 467
 468   // Analyizing the clusters: FRONT ACTIVE ELEMENTS
 469
 470   Int_t nClustersFront = fMFTClusterArrayFront[planeId]->GetEntries();
 471   AliDebug(2, Form("There are %3d clusters in plane %02d FRONT\n", nClustersFront, planeId));
 472
 473   for (Int_t iCluster=0; iCluster<nClustersFront; iCluster++) {
 474
 475     Bool_t isGoodChi2 = kFALSE;
 476
 477     AliMFTCluster *cluster = (AliMFTCluster*) fMFTClusterArrayFront[planeId]->At(iCluster);
 478     Double_t chi2 = (1./(corrFact*corrFact)) * TryOneCluster(currentParamForResearchFront, cluster);     // describes the compatibility between the track and the cluster
 479     if (chi2<chi2cut) isGoodChi2 = kTRUE;
 480
 481     if (isGoodChi2) {
 482       AliDebug(3, Form("accepting cluster: chi2=%f (cut = %f)\n", chi2, chi2cut));
 483       AliMuonForwardTrack *newTrack = new ((*fCandidateTracks)[fCandidateTracks->GetEntriesFast()]) AliMuonForwardTrack(*fCurrentTrack);
 484       if (fCandidateTracks->GetEntriesFast() > fMaxNCandidates) return kDiverged;
 485       newTrack->AddTrackParamAtMFTCluster(currentParamFront, *cluster);    // creating new track param and attaching the cluster
 486       AliDebug(2, Form("After plane %02d: newTrack->GetNMFTClusters() = %d (fCurrentTrack->GetNMFTClusters() = %d)",
 487                        planeId, newTrack->GetNMFTClusters(), fCurrentTrack->GetNMFTClusters()));
 488       newTrack->SetPlaneExists(planeId);
 489     }
 490     else AliDebug(3, Form("discarding cluster: chi2=%f (cut = %f)\n", chi2, chi2cut));
 491
 492   }
 493
 494   // Analyizing the clusters: BACK ACTIVE ELEMENTS
 495
 496   Int_t nClustersBack = fMFTClusterArrayBack[planeId]->GetEntries();
 497   AliDebug(2, Form("There are %3d clusters in plane %02d BACK\n", nClustersBack, planeId));
 498
 499   for (Int_t iCluster=0; iCluster<nClustersBack; iCluster++) {
 500
 501     Bool_t isGoodChi2 = kFALSE;
 502
 503     AliMFTCluster *cluster = (AliMFTCluster*) fMFTClusterArrayBack[planeId]->At(iCluster);
 504     Double_t chi2 = (1./(corrFact*corrFact)) * TryOneCluster(currentParamForResearchBack, cluster);     // describes the compatibility between the track and the cluster
 505     if (chi2<chi2cut) isGoodChi2 = kTRUE;
 506
 507     if (isGoodChi2) {
 508       AliDebug(3,Form("accepting cluster: chi2=%f (cut = %f)\n", chi2, chi2cut));
 509       AliMuonForwardTrack *newTrack = new ((*fCandidateTracks)[fCandidateTracks->GetEntriesFast()]) AliMuonForwardTrack(*fCurrentTrack);
 510       if (fCandidateTracks->GetEntriesFast() > fMaxNCandidates) return kDiverged;
 511       newTrack->AddTrackParamAtMFTCluster(currentParamBack, *cluster);    // creating new track param and attaching the cluster
 512       AliDebug(2, Form("After plane %02d: newTrack->GetNMFTClusters() = %d (fCurrentTrack->GetNMFTClusters() = %d)",
 513                        planeId, newTrack->GetNMFTClusters(), fCurrentTrack->GetNMFTClusters()));
 514       newTrack->SetPlaneExists(planeId);
 515     }
 516     else AliDebug(3,Form("discarding cluster: chi2=%f (cut = %f)\n", chi2, chi2cut));
 517
 518   }
 519
 520   //---------------------------------------------------------------------------------------------
 521
 522   return kConverged;
 523
 524 }
 525
 526 //==========================================================================================================================================
 527
 528 Double_t AliMFTTrackerMU::TryOneCluster(const AliMUONTrackParam &trackParam, AliMFTCluster *cluster) {
 529
 530   // Test the compatibility between the track and the cluster (using trackParam's covariance matrix):
 531   // return the corresponding Chi2
 532   // assume the track parameters are given at the Z of the cluster
 533
 534   // Set differences between trackParam and cluster in the bending and non bending directions
 535   Double_t dX = cluster->GetX() - trackParam.GetNonBendingCoor();
 536   Double_t dY = cluster->GetY() - trackParam.GetBendingCoor();
 537   AliDebug(3,Form("dX = %f, dY = %f\n", dX, dY));
 538
 539   // Calculate errors and covariances
 540   const TMatrixD& kParamCov = trackParam.GetCovariances();
 541   Double_t sigmaX2 = kParamCov(0,0) + cluster->GetErrX2();
 542   Double_t sigmaY2 = kParamCov(2,2) + cluster->GetErrY2();
 543   AliDebug(3, Form("dX2 = %f, dY2 = %f\n", sigmaX2, sigmaY2));
 544   Double_t covXY   = kParamCov(0,2);
 545   Double_t det     = sigmaX2 * sigmaY2 - covXY * covXY;
 546
 547   // Compute chi2
 548   if (det==0.) return 1.e10;
 549   return (dX*dX*sigmaY2 + dY*dY*sigmaX2 - 2.*dX*dY*covXY) / det;
 550
 551 }
 552
 553 //=========================================================================================================================================
 554
 555 Bool_t AliMFTTrackerMU::IsCorrectMatch(AliMFTCluster *cluster, Int_t labelMC) {
 556
 557   Bool_t result = kFALSE;
 558
 559   // check if the cluster belongs to the correct MC track
 560
 561   for (Int_t iTrack=0; iTrack<cluster->GetNMCTracks(); iTrack++) {
 562     if (cluster->GetMCLabel(iTrack)==labelMC) {
 563       result = kTRUE;
 564       break;
 565     }
 566   }
 567
 568   return result;
 569
 570 }
 571
 572 //======================================================================================================================================
 573
 574 void AliMFTTrackerMU::GetVertexFromMC() {
 575
 576   AliRunLoader *fRunLoader = AliRunLoader::Open("galice.root");
 577   if (!fRunLoader) {
 578     AliError("no run loader found in file galice.root");
 579     return;
 580   }
 581
 582   fRunLoader->CdGAFile();
 583   fRunLoader->LoadgAlice();
 584   fRunLoader->LoadHeader();
 585   fRunLoader->GetEvent(gAlice->GetEvNumber());
 586
 587   TArrayF vtx(3);
 588   fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(vtx);
 589
 590   fXExtrapVertex = gRandom->Gaus(vtx[0], AliMFTConstants::fPrimaryVertexResX);
 591   fYExtrapVertex = gRandom->Gaus(vtx[1], AliMFTConstants::fPrimaryVertexResY);
 592   fZExtrapVertex = gRandom->Gaus(vtx[2], AliMFTConstants::fPrimaryVertexResZ);
 593   fXExtrapVertexError = AliMFTConstants::fXVertexTolerance;
 594   fYExtrapVertexError = AliMFTConstants::fYVertexTolerance;
 595   AliInfo(Form("Set ESD vertex from MC (%f +/- %f,  %f +/- %f,  %f)",
 596                fXExtrapVertex,fXExtrapVertexError,fYExtrapVertex,fYExtrapVertexError,fZExtrapVertex));
 597
 598 }
 599
 600 //======================================================================================================================================
 601