changes in the MagF constructor
[u/mrichter/AliRoot.git] / MUON / MUONFakes.C
CommitLineData
66085093 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>
23035434 11#include <TArrayI.h>
66085093 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
43void Prepare(AliMUONRecoParam *&recoParam, Double_t &sigmaCut);
44TTree* GetESDTree(TFile *esdFile);
45Bool_t TrackMatched(AliMUONTrack &track, AliMUONTrack &trackRef, Float_t &fractionOfMatchCluster, Double_t sigmaCut);
46Bool_t IsRecontructible(AliMUONTrack &track, AliMUONRecoParam &recoParam);
47AliMUONTrack* MatchWithTrackRef(AliESDMuonTrack &muonTrack, AliMUONVTrackStore &trackRefStore,
48 Float_t &fractionOfMatchCluster, Bool_t useLabel, Double_t sigmaCut);
49Int_t RemoveConnectedFakes(AliMUONVTrackStore &fakeTrackStore, AliMUONVTrackStore &trackRefStore, AliMUONRecoParam &recoParam,
50 Bool_t useLabel, Double_t sigmaCut, TH1F &hFractionOfConnectedClusters);
51
52//-----------------------------------------------------------------------
53void 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
23035434 104 Int_t nReconstructibleTracks = 0;
105 Int_t nReconstructedTracks = 0;
106 Int_t nEventsWithTrackReconstructedYet = 0;
66085093 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;
23035434 112 Int_t nTotConnectedTracks = 0;
66085093 113 Int_t nTotAdditionalTracks = 0;
23035434 114 Bool_t trackReconstructedYet;
115 TArrayI eventsWithTrackReconstructedYet(10);
116 TArrayI eventsWithFake(10);
117 TArrayI eventsWithAdditionalFake(10);
66085093 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
23035434 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
66085093 141 // loop over ESD tracks
142 Int_t nTrackerTracks = 0;
23035434 143 trackReconstructedYet = kFALSE;
66085093 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++;
23035434 179 if (!IsRecontructible(*matchedTrackRef,*recoParam)) {
180 trackReconstructedYet = kTRUE;
181 nTotTracksReconstructedYet++;
182 }
66085093 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);
23035434 219 nReconstructedTracks += nTrackerTracks;
220 if (trackReconstructedYet) eventsWithTrackReconstructedYet[nEventsWithTrackReconstructedYet++] = iEvent;
66085093 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
23035434 233 eventsWithFake[nEventsWithFake] = iEvent;
66085093 234 nEventsWithFake++;
235 if (nAdditionalTracks > 0) {
23035434 236 eventsWithAdditionalFake[nEventsWithAdditionalFake] = iEvent;
66085093 237 nEventsWithAdditionalFake++;
238 nTotAdditionalTracks += nAdditionalTracks;
239 hNumberOfAdditionalTracks->Fill(nAdditionalTracks);
240 }
241
242 }
243
244 delete fakeTrackStore;
245
246 } // end of loop over events
23035434 247
248 // total number of connected tracks
249 nTotConnectedTracks = hFractionOfConnectedClusters->GetEntries();
66085093 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;
23035434 323 cout << "- Number of reconstructible tracks: " << nReconstructibleTracks << endl;
324 cout << "- Number of reconstructed tracks: " << nReconstructedTracks << endl;
66085093 325 cout << "- Number of matched tracks: " << nTotMatchedTracks << endl;
23035434 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;
66085093 334 cout << "- Number of fake tracks: " << nTotFakeTracks << endl;
23035434 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;
66085093 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//-----------------------------------------------------------------------
365void 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");
4642ac4b 378 AliMagF* field = new AliMagF("Maps","Maps",1.,1., AliMagF::k5kG);
66085093 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//-----------------------------------------------------------------------
407TTree* 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//-----------------------------------------------------------------------
428Bool_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//-----------------------------------------------------------------------
444Bool_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//-----------------------------------------------------------------------
473AliMUONTrack* 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//-----------------------------------------------------------------------
528Int_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