MUON/MUONFakes.C

   1 #if !defined(__CINT__) || defined(__MAKECINT__)
   2 // ROOT includes
   3 #include <TTree.h>
   4 #include <TFile.h>
   5 #include <TH1.h>
   6 #include <TCanvas.h>
   7 #include <Riostream.h>
   8 #include <TROOT.h>
   9 #include <TClonesArray.h>
  10 #include <TGeoGlobalMagField.h>
  11 #include <TArrayI.h>
  12
  13 // STEER includes
  14 #include "AliLog.h"
  15 #include "AliMagF.h"
  16 #include "AliESDEvent.h"
  17 #include "AliESDMuonTrack.h"
  18 #include "AliESDMuonCluster.h"
  19 #include "AliCDBPath.h"
  20 #include "AliCDBEntry.h"
  21 #include "AliCDBManager.h"
  22
  23 // MUON includes
  24 #include "AliMUONTrack.h"
  25 #include "AliMUONVTrackStore.h"
  26 #include "AliMUONTrackParam.h"
  27 #include "AliMUONTrackExtrap.h"
  28 #include "AliMUONESDInterface.h"
  29 #include "AliMUONRecoCheck.h"
  30 #include "AliMUONVCluster.h"
  31 #include "AliMUONRecoParam.h"
  32 #endif
  33
  34 /// \ingroup macros
  35 /// \file MUONFakes.C
  36 ///
  37 /// \author Ph. Pillot, Subatech, March. 2009
  38 ///
  39 /// Macro to study fake tracks by comparing reconstructed tracks with TrackRefs
  40 /// Results are saved in the root file Fakes.root
  41 /// Results are relevent provided that you use the same recoParams as for the reconstruction
  42
  43 void Prepare(AliMUONRecoParam *&recoParam, Double_t &sigmaCut);
  44 TTree* GetESDTree(TFile *esdFile);
  45 Bool_t TrackMatched(AliMUONTrack &track, AliMUONTrack &trackRef, Float_t &fractionOfMatchCluster, Double_t sigmaCut);
  46 Bool_t IsRecontructible(AliMUONTrack &track, AliMUONRecoParam &recoParam);
  47 AliMUONTrack* MatchWithTrackRef(AliESDMuonTrack &muonTrack, AliMUONVTrackStore &trackRefStore,
  48                                 Float_t &fractionOfMatchCluster, Bool_t useLabel, Double_t sigmaCut);
  49 Int_t RemoveConnectedFakes(AliMUONVTrackStore &fakeTrackStore, AliMUONVTrackStore &trackRefStore, AliMUONRecoParam &recoParam,
  50                            Bool_t useLabel, Double_t sigmaCut, TH1F &hFractionOfConnectedClusters);
  51
  52 //-----------------------------------------------------------------------
  53 void MUONFakes(Bool_t useLabel = kFALSE, Int_t FirstEvent = 0, Int_t LastEvent = -1,
  54                const TString esdFileName = "AliESDs.root", const TString SimDir = "./generated/")
  55 {
  56
  57   //Reset ROOT and connect tree file
  58   gROOT->Reset();
  59
  60   // File for histograms and histogram booking
  61   TFile *histoFile = new TFile("Fakes.root", "RECREATE");
  62
  63   TH1F *hNumberOfTracks = new TH1F("hNumberOfTracks","nb of tracks /evt",20,0.,20.);
  64   TH1F *hNumberOfAdditionalTracks = new TH1F("hNumberOfAdditionalTracks","nb of fake - nb of missing track",20,0.,20.);
  65
  66   TH1F *hNumberOfClusters = new TH1F("hNumberOfClusters","nb of clusters /track",20,0.,20.);
  67   TH1F *hNumberOfClustersM = new TH1F("hNumberOfClustersM","nb of clusters /matched track",20,0.,20.);
  68   TH1F *hNumberOfClustersF = new TH1F("hNumberOfClustersF","nb of clusters /fake track",20,0.,20.);
  69   TH1F *hNumberOfClustersMC = new TH1F("hNumberOfClustersMC","nb of clusters /MC track",20,0.,20.);
  70   TH1F *hFractionOfMatchedClusters = new TH1F("hFractionOfMatchedClusters","nb of matched clusters / nb of clusters",110,0.,1.1);
  71   TH1F *hFractionOfConnectedClusters = new TH1F("hFractionOfConnectedClusters","nb of connected clusters / nb of clusters in fake tracks",110,0.,1.1);
  72
  73   TH1F *hChi2PerDof = new TH1F("hChi2PerDof", "track chi2/d.o.f.", 100, 0., 20.);
  74   TH1F *hChi2PerDofM = new TH1F("hChi2PerDofM", "matched track chi2/d.o.f.", 100, 0., 20.);
  75   TH1F *hChi2PerDofF = new TH1F("hChi2PerDofF", "fake track chi2/d.o.f.", 100, 0., 20.);
  76   TH1F *hP = new TH1F("hP", "Muon P distribution (GeV/c)", 100, 0., 200.);
  77   TH1F *hPM = new TH1F("hPM", "matched track P distribution (GeV/c)", 100, 0., 200.);
  78   TH1F *hPF = new TH1F("hPF", "fake track P distribution (GeV/c)", 100, 0., 200.);
  79   TH1F *hPt = new TH1F("hPt", "Muon Pt distribution (GeV/c)", 100, 0., 20.);
  80   TH1F *hPtM = new TH1F("hPtM", "matched track Pt distribution (GeV/c)", 100, 0., 20.);
  81   TH1F *hPtF = new TH1F("hPtF", "fake track Pt distribution (GeV/c)", 100, 0., 20.);
  82   TH1F *hEta = new TH1F("hEta"," Muon pseudo-rapidity distribution",100,-10,0);
  83   TH1F *hEtaM = new TH1F("hEtaM"," matched track pseudo-rapidity distribution",100,-10,0);
  84   TH1F *hEtaF = new TH1F("hEtaF"," fake track pseudo-rapidity distribution",100,-10,0);
  85   TH1F *hPhi = new TH1F("hPhi"," Muon phi distribution",100,-1.,9.);
  86   TH1F *hPhiM = new TH1F("hPhiM"," matched track phi distribution",100,-1.,9.);
  87   TH1F *hPhiF = new TH1F("hPhiF"," fake track phi distribution",100,-1.,9.);
  88
  89   // prepare for analysis
  90   AliMUONRecoParam *recoParam = 0x0;
  91   Double_t sigmaCut = -1;
  92   Prepare(recoParam, sigmaCut);
  93
  94   // link to reconstructed tracks
  95   TFile* esdFile = TFile::Open(esdFileName);
  96   TTree* esdTree = GetESDTree(esdFile);
  97   AliESDEvent* esd = new AliESDEvent();
  98   esd->ReadFromTree(esdTree);
  99
 100   // link to simulated tracks
 101   AliMUONRecoCheck rc(esdFileName, SimDir);
 102
 103   // initialize global counters
 104   Int_t nReconstructibleTracks = 0;
 105   Int_t nReconstructedTracks = 0;
 106   Int_t nEventsWithTrackReconstructedYet = 0;
 107   Int_t nEventsWithFake = 0;
 108   Int_t nEventsWithAdditionalFake = 0;
 109   Int_t nTotMatchedTracks = 0;
 110   Int_t nTotTracksReconstructedYet = 0;
 111   Int_t nTotFakeTracks = 0;
 112   Int_t nTotConnectedTracks = 0;
 113   Int_t nTotAdditionalTracks = 0;
 114   Bool_t trackReconstructedYet;
 115   TArrayI eventsWithTrackReconstructedYet(10);
 116   TArrayI eventsWithFake(10);
 117   TArrayI eventsWithAdditionalFake(10);
 118
 119   // Loop over ESD events
 120   FirstEvent = TMath::Max(0, FirstEvent);
 121   LastEvent = (LastEvent>=0) ? TMath::Min((Int_t)esdTree->GetEntries() - 1, LastEvent) : (Int_t)esdTree->GetEntries() - 1;
 122   for (Int_t iEvent = FirstEvent; iEvent <= LastEvent; iEvent++) {
 123
 124     // get the ESD of current event
 125     esdTree->GetEvent(iEvent);
 126     if (!esd) {
 127       Error("CheckESD", "no ESD object found for event %d", iEvent);
 128       return;
 129     }
 130
 131     // convert TrackRef to MUON tracks
 132     AliMUONVTrackStore* trackRefStore = rc.TrackRefs(iEvent);
 133
 134     // count the number of reconstructible tracks
 135     TIter next(trackRefStore->CreateIterator());
 136     AliMUONTrack* trackRef;
 137     while ( ( trackRef = static_cast<AliMUONTrack*>(next()) ) ) {
 138       if (IsRecontructible(*trackRef,*recoParam)) nReconstructibleTracks++;
 139     }
 140
 141     // loop over ESD tracks
 142     Int_t nTrackerTracks = 0;
 143     trackReconstructedYet = kFALSE;
 144     AliMUONVTrackStore *fakeTrackStore = AliMUONESDInterface::NewTrackStore();
 145     Int_t nTracks = (Int_t)esd->GetNumberOfMuonTracks() ;
 146     for (Int_t iTrack = 0; iTrack <  nTracks;  iTrack++) {
 147
 148       AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
 149
 150       // skip ghosts
 151       if (!muonTrack->ContainTrackerData()) continue;
 152       nTrackerTracks++;
 153
 154       // get track info
 155       Int_t nClusters = muonTrack->GetNClusters();
 156       Double_t normalizedChi2 = muonTrack->GetChi2() / (2. * muonTrack->GetNHit() - 5);
 157       Double_t p = muonTrack->P();
 158       Double_t pT = muonTrack->Pt();
 159       Double_t eta = muonTrack->Eta();
 160       Double_t phi = muonTrack->Phi();
 161
 162       // fill global histograms
 163       hNumberOfClusters->Fill(nClusters);
 164       hChi2PerDof->Fill(normalizedChi2);
 165       hP->Fill(p);
 166       hPt->Fill(pT);
 167       hEta->Fill(eta);
 168       hPhi->Fill(phi);
 169
 170       // try to match the reconstructed track with a simulated one
 171       Float_t fractionOfMatchCluster = 0.;
 172       AliMUONTrack* matchedTrackRef = MatchWithTrackRef(*muonTrack, *trackRefStore, fractionOfMatchCluster, useLabel, sigmaCut);
 173
 174       // take actions according to matching result
 175       if (matchedTrackRef) {
 176
 177         // global counter
 178         nTotMatchedTracks++;
 179         if (!IsRecontructible(*matchedTrackRef,*recoParam)) {
 180           trackReconstructedYet = kTRUE;
 181           nTotTracksReconstructedYet++;
 182         }
 183
 184         // fill histograms
 185         hFractionOfMatchedClusters->Fill(fractionOfMatchCluster);
 186         hNumberOfClustersMC->Fill(matchedTrackRef->GetNClusters());
 187         hNumberOfClustersM->Fill(nClusters);
 188         hChi2PerDofM->Fill(normalizedChi2);
 189         hPM->Fill(p);
 190         hPtM->Fill(pT);
 191         hEtaM->Fill(eta);
 192         hPhiM->Fill(phi);
 193
 194         // remove already matched trackRefs
 195         trackRefStore->Remove(*matchedTrackRef);
 196
 197       } else {
 198
 199         // global counter
 200         nTotFakeTracks++;
 201
 202         // fill histograms
 203         hNumberOfClustersF->Fill(nClusters);
 204         hChi2PerDofF->Fill(normalizedChi2);
 205         hPF->Fill(p);
 206         hPtF->Fill(pT);
 207         hEtaF->Fill(eta);
 208         hPhiF->Fill(phi);
 209
 210         // store fake tracks
 211         AliMUONESDInterface::Add(*muonTrack, *fakeTrackStore);
 212
 213       }
 214
 215     } // end of loop over ESD tracks
 216
 217     // fill histograms
 218     hNumberOfTracks->Fill(nTrackerTracks);
 219     nReconstructedTracks += nTrackerTracks;
 220     if (trackReconstructedYet) eventsWithTrackReconstructedYet[nEventsWithTrackReconstructedYet++] = iEvent;
 221
 222     // count the number the additional fake tracks
 223     if (fakeTrackStore->GetSize() > 0) {
 224
 225       // remove the most connected fake tracks
 226       Int_t nFreeMissingTracks = RemoveConnectedFakes(*fakeTrackStore, *trackRefStore, *recoParam,
 227                                                       useLabel, sigmaCut, *hFractionOfConnectedClusters);
 228
 229       // remove the remaining free reconstructible tracks
 230       Int_t nAdditionalTracks = fakeTrackStore->GetSize() - nFreeMissingTracks;
 231
 232       // fill histograms
 233       eventsWithFake[nEventsWithFake] = iEvent;
 234       nEventsWithFake++;
 235       if (nAdditionalTracks > 0) {
 236         eventsWithAdditionalFake[nEventsWithAdditionalFake] = iEvent;
 237         nEventsWithAdditionalFake++;
 238         nTotAdditionalTracks += nAdditionalTracks;
 239         hNumberOfAdditionalTracks->Fill(nAdditionalTracks);
 240       }
 241
 242     }
 243
 244     delete fakeTrackStore;
 245
 246   } // end of loop over events
 247
 248   // total number of connected tracks
 249   nTotConnectedTracks = hFractionOfConnectedClusters->GetEntries();
 250
 251   // plot results
 252   TCanvas cFakesSummary("cFakesSummary","cFakesSummary",900,600);
 253   cFakesSummary.Divide(3,2);
 254   cFakesSummary.cd(1);
 255   cFakesSummary.GetPad(1)->SetLogy();
 256   hNumberOfClusters->Draw();
 257   hNumberOfClusters->SetMinimum(0.5);
 258   hNumberOfClustersM->Draw("same");
 259   hNumberOfClustersM->SetLineColor(4);
 260   hNumberOfClustersF->Draw("same");
 261   hNumberOfClustersF->SetLineColor(2);
 262   hNumberOfClustersF->SetFillColor(2);
 263   hNumberOfClustersF->SetFillStyle(3017);
 264   cFakesSummary.cd(2);
 265   cFakesSummary.GetPad(2)->SetLogy();
 266   hChi2PerDof->Draw();
 267   hChi2PerDof->SetMinimum(0.5);
 268   hChi2PerDofM->Draw("same");
 269   hChi2PerDofM->SetLineColor(4);
 270   hChi2PerDofF->Draw("same");
 271   hChi2PerDofF->SetLineColor(2);
 272   hChi2PerDofF->SetFillColor(2);
 273   hChi2PerDofF->SetFillStyle(3017);
 274   cFakesSummary.cd(3);
 275   cFakesSummary.GetPad(3)->SetLogy();
 276   hP->Draw();
 277   hP->SetMinimum(0.5);
 278   hPM->Draw("same");
 279   hPM->SetLineColor(4);
 280   hPF->Draw("same");
 281   hPF->SetLineColor(2);
 282   hPF->SetFillColor(2);
 283   hPF->SetFillStyle(3017);
 284   cFakesSummary.cd(4);
 285   cFakesSummary.GetPad(4)->SetLogy();
 286   hPt->Draw();
 287   hPt->SetMinimum(0.5);
 288   hPtM->Draw("same");
 289   hPtM->SetLineColor(4);
 290   hPtF->Draw("same");
 291   hPtF->SetLineColor(2);
 292   hPtF->SetFillColor(2);
 293   hPtF->SetFillStyle(3017);
 294   cFakesSummary.cd(5);
 295   cFakesSummary.GetPad(5)->SetLogy();
 296   hEta->Draw();
 297   hEta->SetMinimum(0.5);
 298   hEtaM->Draw("same");
 299   hEtaM->SetLineColor(4);
 300   hEtaF->Draw("same");
 301   hEtaF->SetLineColor(2);
 302   hEtaF->SetFillColor(2);
 303   hEtaF->SetFillStyle(3017);
 304   cFakesSummary.cd(6);
 305   cFakesSummary.GetPad(6)->SetLogy();
 306   hPhi->Draw();
 307   hPhi->SetMinimum(0.5);
 308   hPhiM->Draw("same");
 309   hPhiM->SetLineColor(4);
 310   hPhiF->Draw("same");
 311   hPhiF->SetLineColor(2);
 312   hPhiF->SetFillColor(2);
 313   hPhiF->SetFillStyle(3017);
 314
 315   // save results
 316   histoFile->cd();
 317   histoFile->Write();
 318   cFakesSummary.Write();
 319   histoFile->Close();
 320
 321   // print results
 322   cout << endl;
 323   cout << "- Number of reconstructible tracks: " << nReconstructibleTracks << endl;
 324   cout << "- Number of reconstructed tracks: " << nReconstructedTracks << endl;
 325   cout << "- Number of matched tracks: " << nTotMatchedTracks << endl;
 326   cout << "  (including " << nTotTracksReconstructedYet << " track(s) matched with a TrackRef that is not reconstructible";
 327   if (nTotTracksReconstructedYet > 0) {
 328     for(Int_t i=0; i<nEventsWithTrackReconstructedYet; i++){
 329       if (i==0) cout << " (eventID = " << eventsWithTrackReconstructedYet[i];
 330       else cout << ", " << eventsWithTrackReconstructedYet[i];
 331     }
 332     cout << "))" << endl;
 333   } else cout << ")" << endl;
 334   cout << "- Number of fake tracks: " << nTotFakeTracks << endl;
 335   cout << "  (including " << nTotConnectedTracks << " track(s) still connected to a reconstructible one)" << endl;
 336   cout << "  (including " << nTotAdditionalTracks << " additional track(s) (compared to the number of expected ones))" << endl;
 337   cout << "- Number of events with fake track(s): " << nEventsWithFake;
 338   if (nEventsWithFake > 0) {
 339     for(Int_t i=0; i<nEventsWithFake; i++){
 340       if (i==0) cout << " (eventID = " << eventsWithFake[i];
 341       else cout << ", " << eventsWithFake[i];
 342     }
 343     cout << ")" << endl;
 344   } else cout << endl;
 345   cout << "  (including " << nEventsWithAdditionalFake << " events with additional track(s)";
 346   if (nEventsWithAdditionalFake > 0) {
 347     for(Int_t i=0; i<nEventsWithAdditionalFake; i++){
 348       if (i==0) cout << " (eventID = " << eventsWithAdditionalFake[i];
 349       else cout << ", " << eventsWithAdditionalFake[i];
 350     }
 351     cout << "))" << endl;
 352   } else cout << ")" << endl;
 353   cout << endl;
 354   cout << "REMINDER: results are relevent provided that you use the same recoParams as for the reconstruction" << endl;
 355   cout << endl;
 356
 357   // clear memory
 358   delete esd;
 359   esdFile->Close();
 360   delete recoParam;
 361
 362 }
 363
 364 //-----------------------------------------------------------------------
 365 void Prepare(AliMUONRecoParam *&recoParam, Double_t &sigmaCut)
 366 {
 367   /// Set the magnetic field and return recoParam and sigmaCut to associate cluster and trackRef
 368
 369   // prepare OCDB access
 370   AliCDBManager* man = AliCDBManager::Instance();
 371   man->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
 372   man->SetRun(0);
 373
 374   // set  mag field
 375   // waiting for mag field in CDB
 376   if (!TGeoGlobalMagField::Instance()->GetField()) {
 377     printf("Loading field map...\n");
 378     AliMagF* field = new AliMagF("Maps","Maps",1.,1., AliMagF::k5kG);
 379     TGeoGlobalMagField::Instance()->SetField(field);
 380   }
 381   // set the magnetic field for track extrapolations
 382   AliMUONTrackExtrap::SetField();
 383
 384   // Load initial reconstruction parameters from OCDB
 385   AliCDBPath path("MUON","Calib","RecoParam");
 386   AliCDBEntry *entry=man->Get(path.GetPath());
 387   if(entry) {
 388     recoParam = dynamic_cast<AliMUONRecoParam*>(entry->GetObject());
 389     entry->SetOwner(0);
 390     AliCDBManager::Instance()->UnloadFromCache(path.GetPath());
 391   }
 392   if (!recoParam) {
 393     printf("Couldn't find RecoParam object in OCDB: create default one");
 394     recoParam = AliMUONRecoParam::GetLowFluxParam();
 395   }
 396
 397   Info("MUONFakes", "\n recontruction parameters:");
 398   recoParam->Print("FULL");
 399   AliMUONESDInterface::ResetTracker(recoParam);
 400
 401   // sigma cut to associate clusters with TrackRefs in case the label are not used
 402   sigmaCut = (recoParam->ImproveTracks()) ? recoParam->GetSigmaCutForImprovement() : recoParam->GetSigmaCutForTracking();
 403
 404 }
 405
 406 //-----------------------------------------------------------------------
 407 TTree* GetESDTree(TFile *esdFile)
 408 {
 409   /// Check that the file is properly open
 410   /// Return pointer to the ESD Tree
 411
 412   if (!esdFile || !esdFile->IsOpen()) {
 413     Error("GetESDTree", "opening ESD file failed");
 414     exit(-1);
 415   }
 416
 417   TTree* tree = (TTree*) esdFile->Get("esdTree");
 418   if (!tree) {
 419     Error("GetESDTree", "no ESD tree found");
 420     exit(-1);
 421   }
 422
 423   return tree;
 424
 425 }
 426
 427 //-----------------------------------------------------------------------
 428 Bool_t TrackMatched(AliMUONTrack &track, AliMUONTrack &trackRef, Float_t &fractionOfMatchCluster, Double_t sigmaCut)
 429 {
 430   /// Try to match 2 tracks
 431
 432   Bool_t compTrack[10];
 433   Int_t nMatchClusters = track.CompatibleTrack(&trackRef, sigmaCut, compTrack);
 434   fractionOfMatchCluster = ((Float_t)nMatchClusters) / ((Float_t)track.GetNClusters());
 435
 436   if ((compTrack[0] || compTrack[1] || compTrack[2] || compTrack[3]) && // at least 1 cluster matched in st 1 & 2
 437       (compTrack[6] || compTrack[7] || compTrack[8] || compTrack[9]) && // at least 1 cluster matched in st 4 & 5
 438       fractionOfMatchCluster > 0.5) return kTRUE;                       // more than 50% of clusters matched
 439   else return kFALSE;
 440
 441 }
 442
 443 //-----------------------------------------------------------------------
 444 Bool_t IsRecontructible(AliMUONTrack &track, AliMUONRecoParam &recoParam)
 445 {
 446   /// Check il the track is reconstructible
 447   Int_t nMinChHitInSt45 = (recoParam.MakeMoreTrackCandidates()) ? 2 : 3;
 448   Int_t currentCh, previousCh = -1, nChHitInSt45 = 0;
 449   Bool_t clusterInSt[5];
 450   for (Int_t iSt = 0; iSt < 5; iSt++) clusterInSt[iSt] = !recoParam.RequestStation(iSt);
 451
 452   AliMUONTrackParam* trackParam = static_cast<AliMUONTrackParam*>(track.GetTrackParamAtCluster()->First());
 453   while (trackParam) {
 454
 455     currentCh = trackParam->GetClusterPtr()->GetChamberId();
 456
 457     clusterInSt[currentCh/2] = kTRUE;
 458
 459     if (currentCh > 5 && currentCh != previousCh) {
 460       nChHitInSt45++;
 461       previousCh = currentCh;
 462     }
 463
 464     trackParam = static_cast<AliMUONTrackParam*>(track.GetTrackParamAtCluster()->After(trackParam));
 465   }
 466
 467   return (clusterInSt[0] && clusterInSt[1] && clusterInSt[2] &&
 468           clusterInSt[3] && clusterInSt[4] && nChHitInSt45 >= nMinChHitInSt45);
 469
 470 }
 471
 472 //-----------------------------------------------------------------------
 473 AliMUONTrack* MatchWithTrackRef(AliESDMuonTrack &muonTrack, AliMUONVTrackStore &trackRefStore,
 474                                 Float_t &fractionOfMatchCluster, Bool_t useLabel, Double_t sigmaCut)
 475 {
 476   /// Return if the trackRef matched with the reconstructed track and the fraction of matched clusters
 477
 478   AliMUONTrack *matchedTrackRef = 0x0;
 479   fractionOfMatchCluster = 0.;
 480
 481   if (useLabel) { // by using the MC label
 482
 483     // get the corresponding simulated track if any
 484     Int_t label = muonTrack.GetLabel();
 485     matchedTrackRef = (AliMUONTrack*) trackRefStore.FindObject(label);
 486
 487     // get the fraction of matched clusters
 488     if (matchedTrackRef) {
 489       Int_t nMatchClusters = 0;
 490       if (muonTrack.ClustersStored()) {
 491         AliESDMuonCluster* cluster = (AliESDMuonCluster*) muonTrack.GetClusters().First();
 492         while (cluster) {
 493           if (cluster->GetLabel() == label) nMatchClusters++;
 494           cluster = (AliESDMuonCluster*) muonTrack.GetClusters().After(cluster);
 495         }
 496       }
 497       fractionOfMatchCluster = ((Float_t)nMatchClusters) / ((Float_t)muonTrack.GetNClusters());
 498     }
 499
 500   } else { // by comparing cluster/TrackRef positions
 501
 502     // convert ESD track to MUON track
 503     AliMUONTrack track;
 504     AliMUONESDInterface::ESDToMUON(muonTrack,track);
 505
 506     // look for the corresponding simulated track if any
 507     TIter next(trackRefStore.CreateIterator());
 508     AliMUONTrack* trackRef;
 509     while ( ( trackRef = static_cast<AliMUONTrack*>(next()) ) ) {
 510
 511       // check compatibility
 512       Float_t f = 0.;
 513       if (TrackMatched(track, *trackRef, f, sigmaCut)) {
 514         matchedTrackRef = trackRef;
 515         fractionOfMatchCluster = f;
 516         break;
 517       }
 518
 519     }
 520
 521   }
 522
 523   return matchedTrackRef;
 524
 525 }
 526
 527 //-----------------------------------------------------------------------
 528 Int_t RemoveConnectedFakes(AliMUONVTrackStore &fakeTrackStore, AliMUONVTrackStore &trackRefStore, AliMUONRecoParam &recoParam,
 529                           Bool_t useLabel, Double_t sigmaCut, TH1F &hFractionOfConnectedClusters)
 530 {
 531   /// loop over reconstructible TrackRef not associated with reconstructed track:
 532   /// for each of them, find and remove the most connected the fake track, if any,
 533   /// and fill the histograms with the fraction of connected clusters.
 534   /// Return the number of reconstructible track not connected to any fake
 535
 536   Int_t nFreeMissingTracks = 0;
 537
 538   // loop over trackRefs
 539   TIter next(trackRefStore.CreateIterator());
 540   AliMUONTrack* trackRef;
 541   while ( ( trackRef = static_cast<AliMUONTrack*>(next()) ) ) {
 542
 543     // skip not reconstructible trackRefs
 544     if (!IsRecontructible(*trackRef,recoParam)) continue;
 545
 546     Int_t label = trackRef->GetUniqueID();
 547
 548     // look for the most connected fake track
 549     AliMUONTrack *connectedFake = 0x0;
 550     Float_t fractionOfConnectedClusters = 0.;
 551     TIter next2(fakeTrackStore.CreateIterator());
 552     AliMUONTrack* fakeTrack;
 553     while ( ( fakeTrack = static_cast<AliMUONTrack*>(next2()) ) ) {
 554
 555       // get the number of connected clusters
 556       Int_t nConnectedClusters = 0;
 557       if (useLabel) { // by using the MC label
 558         for (Int_t iCl = 0; iCl < fakeTrack->GetNClusters(); iCl++)
 559           if (((AliMUONTrackParam*) fakeTrack->GetTrackParamAtCluster()->UncheckedAt(iCl))->GetClusterPtr()->GetMCLabel() == label)
 560             nConnectedClusters++;
 561       } else { // by comparing cluster/TrackRef positions
 562         Bool_t compTrack[10];
 563         nConnectedClusters = fakeTrack->CompatibleTrack(trackRef, sigmaCut, compTrack);
 564       }
 565
 566       // skip non-connected fake tracks
 567       if (nConnectedClusters == 0) continue;
 568
 569       // check if it is the most connected fake track
 570       Float_t f = ((Float_t)nConnectedClusters) / ((Float_t)fakeTrack->GetNClusters());
 571       if (f > fractionOfConnectedClusters) {
 572         connectedFake = fakeTrack;
 573         fractionOfConnectedClusters = f;
 574       }
 575
 576     }
 577
 578     // remove the most connected fake track
 579     if (connectedFake) {
 580       hFractionOfConnectedClusters.Fill(fractionOfConnectedClusters);
 581       fakeTrackStore.Remove(*connectedFake);
 582     } else nFreeMissingTracks++;
 583
 584   }
 585
 586   return nFreeMissingTracks;
 587
 588 }
 589