CbmRoot
CbmBinnedTrackerQANew.cxx
Go to the documentation of this file.
1 /*
2  * To change this license header, choose License Headers in Project Properties.
3  * To change this template file, choose Tools | Templates
4  * and open the template in the editor.
5  */
6 
8 #include "CbmMCDataManager.h"
9 #include "CbmMCTrack.h"
10 #include "CbmMatch.h"
11 #include "CbmMuchCluster.h"
12 #include "CbmMuchDigi.h"
13 #include "CbmMuchGeoScheme.h"
14 #include "CbmMuchPixelHit.h"
15 #include "CbmMuchPoint.h"
16 #include "CbmStsCluster.h"
17 #include "CbmStsDigi.h"
18 #include "CbmStsHit.h"
19 #include "CbmStsPoint.h"
20 #include "CbmStsSetup.h"
21 #include "CbmTofPoint.h"
22 #include "CbmTrdCluster.h"
23 #include "FairRun.h"
24 #include "FairRuntimeDb.h"
25 #include "global/CbmGlobalTrack.h"
26 
27 using namespace std;
28 
29 struct MCTrackDesc {
30  map<pair<ECbmModuleId, Int_t>, pair<UInt_t, set<Int_t>>>
31  mcPointMap; // = map<pair<ECbmModuleId, Int_t>, pair<UInt_t, set<Int_t> > > ();// Map an MC point index to TRD station number and a set of reconstructed global track indices.
32  Int_t pdg = -1;
33  Int_t parentInd = -1;
34 };
35 
36 static vector<vector<MCTrackDesc>> gMCTracks;
37 
39 
41  CbmStsSetup* stsSetup = CbmStsSetup::Instance();
42 
43  if (!stsSetup->IsInit()) stsSetup->Init();
44 
45  fSettings = CbmBinnedSettings::Instance();
46 
47  fMinTrackLength +=
48  fSettings->Use(ECbmModuleId::kSts) ? fSettings->GetNofStsStations() : 0;
49  fMinTrackLength +=
50  fSettings->Use(ECbmModuleId::kMuch) ? fSettings->GetNofMuchStations() : 0;
51  fMinTrackLength +=
52  fSettings->Use(ECbmModuleId::kTrd) ? fSettings->GetNofTrdStations() : 0;
53  fMinTrackLength += fSettings->Use(ECbmModuleId::kTof) ? 1 : 0;
54  int canSkipHits = 0.3 * fMinTrackLength;
55  fMinTrackLength -= canSkipHits;
56 
57  FairRootManager* ioman = FairRootManager::Instance();
58 
59  if (0 == ioman) LOG(fatal) << "No FairRootManager";
60 
61  CbmMCDataManager* mcManager =
62  static_cast<CbmMCDataManager*>(ioman->GetObject("MCDataManager"));
63 
64  if (0 == mcManager) LOG(fatal) << "No MC data manager";
65 
66  fMCTracks = mcManager->InitBranch("MCTrack");
67 
68  if (0 == fMCTracks) LOG(fatal) << "No MC tracks in the input file";
69 
70  for (int i = 0; fMCTracks->Size(0, i) >= 0; ++i) {
71  gMCTracks.push_back(vector<MCTrackDesc>());
72 
73  auto nofMcTracks = fMCTracks->Size(0, i);
74  auto& eventTracks = gMCTracks.back();
75  eventTracks.resize(nofMcTracks);
76 
77  for (Int_t j = 0; j < nofMcTracks; ++j) {
78  MCTrackDesc& track = eventTracks[j];
79  const CbmMCTrack* mcTrack =
80  static_cast<const CbmMCTrack*>(fMCTracks->Get(0, i, j));
81  track.pdg = mcTrack->GetPdgCode();
82  track.parentInd = mcTrack->GetMotherId();
83  }
84  } // for (Int_t i = 0; fMCTracks->Size(0, i) >= 0; ++i)
85 
86  fGlobalTracks = static_cast<TClonesArray*>(ioman->GetObject("GlobalTrack"));
87 
88  if (0 == fGlobalTracks) LOG(fatal) << "No global tracks in the input file";
89 
90  if (fSettings->Use(ECbmModuleId::kSts)) {
91  fStsTracks = static_cast<TClonesArray*>(ioman->GetObject("StsTrack"));
92 
93  if (0 == fStsTracks) LOG(fatal) << "No sts tracks in the input file";
94 
95  fStsHits = static_cast<TClonesArray*>(ioman->GetObject("StsHit"));
96 
97  if (0 == fStsHits) LOG(fatal) << "No sts hits in the input file";
98 
99  fStsClusters = static_cast<TClonesArray*>(ioman->GetObject("StsCluster"));
100 
101  if (0 == fStsClusters) LOG(fatal) << "No sts clusters in the input file";
102 
103  fStsDigis = static_cast<TClonesArray*>(ioman->GetObject("StsDigi"));
104 
105  if (0 == fStsDigis) LOG(fatal) << "No sts digis in the input file";
106 
107  fStsDigiMatches =
108  static_cast<TClonesArray*>(ioman->GetObject("StsDigiMatch"));
109 
110  if (0 == fStsDigiMatches)
111  LOG(fatal) << "No sts digi matches in the input file";
112 
113  fStsPoints = mcManager->InitBranch("StsPoint");
114 
115  if (0 == fStsPoints) LOG(fatal) << "No sts MC points in the input file";
116 
117  for (Int_t i = 0; fStsPoints->Size(0, i) >= 0; ++i) {
118  Int_t nofPoints = fStsPoints->Size(0, i);
119  vector<MCTrackDesc>& tracks = gMCTracks[i];
120 
121  for (Int_t j = 0; j < nofPoints; ++j) {
122  const CbmStsPoint* stsPoint =
123  static_cast<const CbmStsPoint*>(fStsPoints->Get(0, i, j));
124  Int_t trackId = stsPoint->GetTrackID();
125  Int_t stationNumber =
126  CbmStsSetup::Instance()->GetStationNumber(stsPoint->GetDetectorID());
127  auto& trackDesk = tracks[trackId];
128  trackDesk.mcPointMap[make_pair(ECbmModuleId::kSts, j)] = {stationNumber,
129  set<Int_t>()};
130  }
131  }
132  } // if (fSettings->Use(kSts))
133 
134  if (fSettings->Use(ECbmModuleId::kMuch)) {
135  fMuchTracks = static_cast<TClonesArray*>(ioman->GetObject("MuchTrack"));
136 
137  if (0 == fMuchTracks) LOG(fatal) << "No much tracks in the input file";
138 
139  fMuchHits = static_cast<TClonesArray*>(ioman->GetObject("MuchPixelHit"));
140 
141  if (0 == fMuchHits) LOG(fatal) << "No much hits in the input file";
142 
143  fMuchClusters = static_cast<TClonesArray*>(ioman->GetObject("MuchCluster"));
144 
145  if (0 == fMuchClusters) LOG(fatal) << "No much clusters in the input file";
146 
147  fMuchDigis = static_cast<TClonesArray*>(ioman->GetObject("MuchDigi"));
148 
149  if (0 == fMuchDigis) LOG(fatal) << "No much digis in the input file";
150 
151  fMuchDigiMatches =
152  static_cast<TClonesArray*>(ioman->GetObject("MuchDigiMatch"));
153 
154  if (0 == fMuchDigiMatches)
155  LOG(fatal) << "No much digi matches in the input file";
156 
157  fMuchPoints = mcManager->InitBranch("MuchPoint");
158 
159  if (0 == fMuchPoints) LOG(fatal) << "No much MC points in the input file";
160 
161  for (Int_t i = 0; fMuchPoints->Size(0, i) >= 0; ++i) {
162  Int_t nofPoints = fMuchPoints->Size(0, i);
163  vector<MCTrackDesc>& tracks = gMCTracks[i];
164 
165  for (Int_t j = 0; j < nofPoints; ++j) {
166  const CbmMuchPoint* muchPoint =
167  static_cast<const CbmMuchPoint*>(fMuchPoints->Get(0, i, j));
168  Int_t trackId = muchPoint->GetTrackID();
169  int muchStationNumber =
170  CbmMuchGeoScheme::GetStationIndex(muchPoint->GetDetectorID());
171  int layerNumber =
172  CbmMuchGeoScheme::GetLayerIndex(muchPoint->GetDetectorID());
173  int stationNumber = muchStationNumber * 3 + layerNumber;
174  auto& trackDesk = tracks[trackId];
175  trackDesk.mcPointMap[make_pair(ECbmModuleId::kMuch, j)] = {
176  stationNumber, set<Int_t>()};
177  }
178  }
179  } // if (fSettings->Use(kMuch))
180 
181  if (fSettings->Use(ECbmModuleId::kTrd)) {
182  fTrdTracks = static_cast<TClonesArray*>(ioman->GetObject("TrdTrack"));
183 
184  if (0 == fTrdTracks) LOG(fatal) << "No trd tracks in the input file";
185 
186  fTrdHits = static_cast<TClonesArray*>(ioman->GetObject("TrdHit"));
187 
188  if (0 == fTrdHits) LOG(fatal) << "No trd hits in the input file";
189 
190  fTrdClusters = static_cast<TClonesArray*>(ioman->GetObject("TrdCluster"));
191 
192  if (0 == fTrdClusters) LOG(fatal) << "No global tracks in the input file";
193 
194  fTrdDigiMatches =
195  static_cast<TClonesArray*>(ioman->GetObject("TrdDigiMatch"));
196 
197  if (0 == fTrdDigiMatches)
198  LOG(fatal) << "No trd hit to digi matches in the input file";
199 
200  fTrdPoints = mcManager->InitBranch("TrdPoint");
201 
202  if (0 == fTrdPoints) LOG(fatal) << "No trd MC points in the input file";
203 
204  for (Int_t i = 0; fTrdPoints->Size(0, i) >= 0; ++i) {
205  Int_t nofPoints = fTrdPoints->Size(0, i);
206  vector<MCTrackDesc>& tracks = gMCTracks[i];
207 
208  for (Int_t j = 0; j < nofPoints; ++j) {
209  const CbmTrdPoint* trdPoint =
210  static_cast<const CbmTrdPoint*>(fTrdPoints->Get(0, i, j));
211  Int_t trackId = trdPoint->GetTrackID();
212  int stationNumber =
214  auto& trackDesk = tracks[trackId];
215  trackDesk.mcPointMap[make_pair(ECbmModuleId::kTrd, j)] = {stationNumber,
216  set<Int_t>()};
217  }
218  }
219  } // if (fSettings->Use(kTrd))
220 
221  if (fSettings->Use(ECbmModuleId::kTof)) {
222  fTofHits = static_cast<TClonesArray*>(ioman->GetObject("TofHit"));
223 
224  if (0 == fTofHits) LOG(fatal) << "No tof hits in the input file";
225 
226  fTofHitDigiMatches =
227  static_cast<TClonesArray*>(ioman->GetObject("TofDigiMatch"));
228 
229  if (0 == fTofHitDigiMatches)
230  LOG(fatal) << "No tof hit to digi matches in the input file";
231 
232  fTofDigiPointMatches =
233  static_cast<TClonesArray*>(ioman->GetObject("TofDigiMatchPoints"));
234 
235  if (0 == fTofDigiPointMatches)
236  LOG(fatal) << "No tof digi to point matches in the input file";
237 
238  fTofDigis = static_cast<TClonesArray*>(ioman->GetObject("TofDigi"));
239 
240  if (0 == fTofDigis) LOG(fatal) << "No tof digis in the input file";
241 
242  fTofPoints = mcManager->InitBranch("TofPoint");
243 
244  if (0 == fTofPoints) LOG(fatal) << "No tof MC points in the input file";
245 
246  for (Int_t i = 0; fTofPoints->Size(0, i) >= 0; ++i) {
247  Int_t nofPoints = fTofPoints->Size(0, i);
248  vector<MCTrackDesc>& tracks = gMCTracks[i];
249 
250  for (Int_t j = 0; j < nofPoints; ++j) {
251  const CbmTofPoint* tofPoint =
252  static_cast<const CbmTofPoint*>(fTofPoints->Get(0, i, j));
253  Int_t trackId = tofPoint->GetTrackID();
254  auto& trackDesk = tracks[trackId];
255  trackDesk.mcPointMap[make_pair(ECbmModuleId::kTof, j)] = {0,
256  set<Int_t>()};
257  }
258  }
259  } // if (fSettings->Use(kTof))
260 
261  return kSUCCESS;
262 }
263 
264 static Int_t gEventNumber = 0;
265 static int gNofRecoTracks = 0;
266 static int gNofMatchedRecoTracks = 0;
267 static int gNofClones = 0;
268 
270  auto& eventMCTracks = gMCTracks[gEventNumber];
271  Int_t nofGlobalTracks = fGlobalTracks->GetEntriesFast();
272 
273  for (Int_t i = 0; i < nofGlobalTracks; ++i) {
274  const CbmGlobalTrack* globalTrack =
275  static_cast<const CbmGlobalTrack*>(fGlobalTracks->At(i));
276  Int_t stsTrackIndex = globalTrack->GetStsTrackIndex();
277  Int_t muchTrackIndex = globalTrack->GetMuchTrackIndex();
278  Int_t trdTrackIndex = globalTrack->GetTrdTrackIndex();
279  Int_t tofHitIndex = globalTrack->GetTofHitIndex();
280 
281  map<void*, set<pair<ECbmModuleId, Int_t>>> mcTrackPtrs;
282  Int_t nofTrackHits = 0;
283 
284  if (fSettings->Use(ECbmModuleId::kSts) && stsTrackIndex >= 0) {
285  const CbmStsTrack* stsTrack =
286  static_cast<const CbmStsTrack*>(fStsTracks->At(stsTrackIndex));
287  Int_t nofStsHits = stsTrack->GetNofHits();
288  nofTrackHits += nofStsHits;
289 
290  for (Int_t j = 0; j < nofStsHits; ++j) {
291  Int_t stsHitInd = stsTrack->GetStsHitIndex(j);
292  const CbmStsHit* stsHit =
293  static_cast<const CbmStsHit*>(fStsHits->At(stsHitInd));
294  //Int_t stationNumber = CbmStsSetup::Instance()->GetStationNumber(stsHit->GetAddress());
295  Int_t frontClusterInd = stsHit->GetFrontClusterId();
296  Int_t backClusterInd = stsHit->GetBackClusterId();
297  const CbmStsCluster* frontCluster =
298  static_cast<const CbmStsCluster*>(fStsClusters->At(frontClusterInd));
299  Int_t nofFrontDigis = frontCluster->GetNofDigis();
300 
301  for (Int_t k = 0; k < nofFrontDigis; ++k) {
302  Int_t stsDigiInd = frontCluster->GetDigi(k);
303  //const CbmStsDigi* stsDigi = static_cast<const CbmStsDigi*> (fStsDigis->At(stsDigiInd));
304  const CbmMatch* stsDigiMatch =
305  static_cast<const CbmMatch*>(fStsDigiMatches->At(stsDigiInd));
306  Int_t nofLinks = stsDigiMatch->GetNofLinks();
307 
308  for (Int_t m = 0; m < nofLinks; ++m) {
309  const CbmLink& link = stsDigiMatch->GetLink(m);
310  Int_t eventId = link.GetEntry();
311  Int_t mcPointId = link.GetIndex();
312  const CbmStsPoint* stsPoint = static_cast<const CbmStsPoint*>(
313  fStsPoints->Get(0, eventId, mcPointId));
314  Int_t trackId = stsPoint->GetTrackID();
315  auto& mcTrack = eventMCTracks[trackId];
316  auto recoIter =
317  mcTrack.mcPointMap.find(make_pair(ECbmModuleId::kSts, mcPointId));
318 
319  if (recoIter != mcTrack.mcPointMap.end())
320  recoIter->second.second.insert(i);
321 
322  mcTrackPtrs[&mcTrack].insert(make_pair(ECbmModuleId::kSts, j));
323  }
324  }
325 
326  const CbmStsCluster* backCluster =
327  static_cast<const CbmStsCluster*>(fStsClusters->At(backClusterInd));
328  Int_t nofBackDigis = backCluster->GetNofDigis();
329 
330  for (Int_t k = 0; k < nofBackDigis; ++k) {
331  Int_t stsDigiInd = backCluster->GetDigi(k);
332  //const CbmStsDigi* stsDigi = static_cast<const CbmStsDigi*> (fStsDigis->At(stsDigiInd));
333  const CbmMatch* stsDigiMatch =
334  static_cast<const CbmMatch*>(fStsDigiMatches->At(stsDigiInd));
335  Int_t nofLinks = stsDigiMatch->GetNofLinks();
336 
337  for (Int_t m = 0; m < nofLinks; ++m) {
338  const CbmLink& link = stsDigiMatch->GetLink(m);
339  Int_t eventId = link.GetEntry();
340  Int_t mcPointId = link.GetIndex();
341  const CbmStsPoint* stsPoint = static_cast<const CbmStsPoint*>(
342  fStsPoints->Get(0, eventId, mcPointId));
343  Int_t trackId = stsPoint->GetTrackID();
344  auto& mcTrack = eventMCTracks[trackId];
345  auto recoIter =
346  mcTrack.mcPointMap.find(make_pair(ECbmModuleId::kSts, mcPointId));
347 
348  if (recoIter != mcTrack.mcPointMap.end())
349  recoIter->second.second.insert(i);
350 
351  mcTrackPtrs[&mcTrack].insert(make_pair(ECbmModuleId::kSts, j));
352  }
353  }
354  } // for (Int_t i = 0; i < nofStsHits; ++i)
355  } // if (fSettings->Use(kSts) && stsTrackIndex >= 0)
356 
357  if (fSettings->Use(ECbmModuleId::kMuch) && muchTrackIndex >= 0) {
358  const CbmMuchTrack* muchTrack =
359  static_cast<const CbmMuchTrack*>(fMuchTracks->At(muchTrackIndex));
360  Int_t nofMuchHits = muchTrack->GetNofHits();
361  nofTrackHits += nofMuchHits;
362 
363  for (Int_t j = 0; j < nofMuchHits; ++j) {
364  Int_t muchHitInd = muchTrack->GetHitIndex(j);
365  const CbmMuchPixelHit* muchHit =
366  static_cast<const CbmMuchPixelHit*>(fMuchHits->At(muchHitInd));
367  //Int_t muchStationNumber = CbmMuchGeoScheme::GetStationIndex(muchHit->GetAddress());
368  //Int_t layerNumber = CbmMuchGeoScheme::GetLayerIndex(muchHit->GetAddress());
369  //Int_t stationNumber = muchStationNumber * 3 + layerNumber;
370  Int_t clusterId = muchHit->GetRefId();
371  const CbmMuchCluster* cluster =
372  static_cast<const CbmMuchCluster*>(fMuchClusters->At(clusterId));
373  Int_t nofDigis = cluster->GetNofDigis();
374 
375  for (Int_t k = 0; k < nofDigis; ++k) {
376  Int_t digiId = cluster->GetDigi(k);
377  //const CbmMuchDigi* digi = static_cast<const CbmMuchDigi*> (fMuchDigis->At(digiId));
378  const CbmMatch* digiMatch =
379  static_cast<const CbmMatch*>(fMuchDigiMatches->At(digiId));
380  Int_t nofLinks = digiMatch->GetNofLinks();
381 
382  for (Int_t m = 0; m < nofLinks; ++m) {
383  const CbmLink& link = digiMatch->GetLink(m);
384  Int_t eventId = link.GetEntry();
385  Int_t mcPointId = link.GetIndex();
386  const CbmMuchPoint* muchPoint = static_cast<const CbmMuchPoint*>(
387  fMuchPoints->Get(0, eventId, mcPointId));
388  Int_t trackId = muchPoint->GetTrackID();
389  auto& mcTrack = eventMCTracks[trackId];
390  auto recoIter = mcTrack.mcPointMap.find(
391  make_pair(ECbmModuleId::kMuch, mcPointId));
392 
393  if (recoIter != mcTrack.mcPointMap.end())
394  recoIter->second.second.insert(i);
395 
396  mcTrackPtrs[&mcTrack].insert(make_pair(ECbmModuleId::kMuch, j));
397  }
398  }
399  }
400  } // if (fSettings->Use(ECbmModuleId::kMuch) && muchTrackIndex >= 0)
401 
402  if (fSettings->Use(ECbmModuleId::kTrd) && trdTrackIndex >= 0) {
403  const CbmTrdTrack* trdTrack =
404  static_cast<const CbmTrdTrack*>(fTrdTracks->At(trdTrackIndex));
405  Int_t nofTrdHits = trdTrack->GetNofHits();
406  nofTrackHits += nofTrdHits;
407 
408  for (Int_t j = 0; j < nofTrdHits; ++j) {
409  Int_t trdHitInd = trdTrack->GetHitIndex(j);
410  const CbmTrdHit* trdHit =
411  static_cast<const CbmTrdHit*>(fTrdHits->At(trdHitInd));
412  Int_t clusterId = trdHit->GetRefId();
413  const CbmTrdCluster* cluster =
414  static_cast<const CbmTrdCluster*>(fTrdClusters->At(clusterId));
415  Int_t nofDigis = cluster->GetNofDigis();
416 
417  for (Int_t k = 0; k < nofDigis; ++k) {
418  Int_t digiId = cluster->GetDigi(k);
419  const CbmMatch* match =
420  static_cast<const CbmMatch*>(fTrdDigiMatches->At(digiId));
421  Int_t nofLinks = match->GetNofLinks();
422 
423  for (Int_t m = 0; m < nofLinks; ++m) {
424  const CbmLink& link = match->GetLink(m);
425  Int_t eventId = link.GetEntry();
426  Int_t mcPointId = link.GetIndex();
427  const CbmTrdPoint* trdPoint = static_cast<const CbmTrdPoint*>(
428  fTrdPoints->Get(0, eventId, mcPointId));
429  Int_t trackId = trdPoint->GetTrackID();
430  auto& mcTrack = eventMCTracks[trackId];
431  auto recoIter =
432  mcTrack.mcPointMap.find(make_pair(ECbmModuleId::kTrd, mcPointId));
433 
434  if (recoIter != mcTrack.mcPointMap.end())
435  recoIter->second.second.insert(i);
436 
437  mcTrackPtrs[&mcTrack].insert(make_pair(ECbmModuleId::kTrd, j));
438  }
439  }
440  } // for (Int_t j = 0; j < nofTrdHits; ++j)
441  } // if (fSettings->Use(kTrd) && trdTrackIndex >= 0)
442 
443  if (fSettings->Use(ECbmModuleId::kTof) && tofHitIndex >= 0) {
444  ++nofTrackHits;
445  const CbmMatch* tofHitMatch =
446  static_cast<const CbmMatch*>(fTofHitDigiMatches->At(tofHitIndex));
447  Int_t nofTofDigis = tofHitMatch->GetNofLinks();
448 
449  for (Int_t j = 0; j < nofTofDigis; ++j) {
450  const CbmLink& digiLink = tofHitMatch->GetLink(j);
451  Int_t digiInd = digiLink.GetIndex();
452  const CbmMatch* pointMatch =
453  static_cast<const CbmMatch*>(fTofDigiPointMatches->At(digiInd));
454  Int_t nofPoints = pointMatch->GetNofLinks();
455 
456  for (Int_t k = 0; k < nofPoints; ++k) {
457  const CbmLink& pointLink = pointMatch->GetLink(k);
458  Int_t mcEventId = pointLink.GetEntry();
459  Int_t mcPointId = pointLink.GetIndex();
460  const CbmTofPoint* tofPoint = static_cast<const CbmTofPoint*>(
461  fTofPoints->Get(0, mcEventId, mcPointId));
462  Int_t trackId = tofPoint->GetTrackID();
463  auto& mcTrack = eventMCTracks[trackId];
464  auto recoIter =
465  mcTrack.mcPointMap.find(make_pair(ECbmModuleId::kTof, mcPointId));
466 
467  if (recoIter != mcTrack.mcPointMap.end())
468  recoIter->second.second.insert(i);
469 
470  mcTrackPtrs[&mcTrack].insert(make_pair(ECbmModuleId::kTof, 0));
471  }
472  }
473  } // if (fSettings->Use(kTof) && tofHitIndex >= 0)
474 
475  ++gNofRecoTracks;
476 
477  for (const auto& j : mcTrackPtrs) {
478  if (j.second.size() >= size_t(0.7 * nofTrackHits)) {
480  break;
481  }
482  }
483  } // for (Int_t i = 0; i < nofGlobalTracks; ++i)
484 
485  ++gEventNumber;
486 } //Exec()
487 
489  int nofRefMCTRacks = 0;
490  int nofMatchedRefMCTracks = 0;
491  uint maxTrackLength = 0;
492 
493  for (const auto& i : gMCTracks) {
494  for (const auto& j : i) {
495  set<pair<ECbmModuleId, UInt_t>> stationNofs;
496 
497  for (const auto& k : j.mcPointMap)
498  stationNofs.insert(make_pair(k.first.first, k.second.first));
499 
500  if (stationNofs.size() > static_cast<size_t>(maxTrackLength))
501  maxTrackLength = stationNofs.size();
502 
503  if (stationNofs.size() < static_cast<size_t>(fMinTrackLength)) continue;
504 
505  map<Int_t, set<pair<ECbmModuleId, UInt_t>>> recoToStationMatches;
506 
507  for (const auto& k : j.mcPointMap) {
508  for (const auto& m : k.second.second)
509  recoToStationMatches[m].insert(
510  make_pair(k.first.first, k.second.first));
511  }
512 
513  int nofMatchedReco = 0;
514 
515  for (const auto& k : recoToStationMatches) {
516  if (k.second.size() >= size_t(stationNofs.size() * 0.7))
517  ++nofMatchedReco;
518  }
519 
520  if (nofMatchedReco > 1) gNofClones += nofMatchedReco - 1;
521 
522  //if (abs(j.pdg) != 11)
523  //continue;
524 
525  if (fPrimaryParticlePdg >= 0
526  && (j.parentInd < 0 || i[j.parentInd].pdg != fPrimaryParticlePdg))
527  continue;
528 
529  ++nofRefMCTRacks;
530 
531  if (nofMatchedReco > 0) ++nofMatchedRefMCTracks;
532  }
533  }
534 
535  cout << "Maximum MC track length: " << maxTrackLength << endl;
536 
537  double res = 100 * nofMatchedRefMCTracks;
538  res /= nofRefMCTRacks;
539  cout << "The reconstruction efficiency: " << res << " % = 100 * "
540  << nofMatchedRefMCTracks << " / " << nofRefMCTRacks << endl;
541 
542  res = 100 * (gNofRecoTracks - gNofMatchedRecoTracks);
543  res /= gNofRecoTracks;
544  cout << "The share of ghosts: " << res << " % = 100 * (" << gNofRecoTracks
545  << " - " << gNofMatchedRecoTracks << ") / " << gNofRecoTracks << endl;
546 
547  res = 100 * gNofClones;
548  res /= gNofRecoTracks;
549  cout << "The share of clones: " << res << " % = 100 * " << gNofClones << " / "
550  << gNofRecoTracks << endl;
551 }
552 
554  fSettings = static_cast<CbmBinnedSettings*>(
555  FairRun::Instance()->GetRuntimeDb()->getContainer("CbmBinnedSettings"));
556 }
557 
CbmMCTrack::GetMotherId
Int_t GetMotherId() const
Definition: CbmMCTrack.h:71
CbmBinnedTrackerQANew::SetParContainers
void SetParContainers()
Definition: CbmBinnedTrackerQANew.cxx:553
CbmBinnedTrackerQANew::Exec
void Exec(Option_t *)
Definition: CbmBinnedTrackerQANew.cxx:269
CbmMatch
Definition: CbmMatch.h:22
CbmMCDataManager.h
CbmMuchDigi.h
CbmMuchPoint
Definition: CbmMuchPoint.h:21
CbmGlobalTrack::GetMuchTrackIndex
Int_t GetMuchTrackIndex() const
Definition: CbmGlobalTrack.h:40
CbmMatch::GetLink
const CbmLink & GetLink(Int_t i) const
Definition: CbmMatch.h:35
CbmMatch::GetNofLinks
Int_t GetNofLinks() const
Definition: CbmMatch.h:38
CbmStsCluster
Data class for STS clusters.
Definition: CbmStsCluster.h:31
CbmTrack::GetNofHits
virtual Int_t GetNofHits() const
Definition: CbmTrack.h:53
CbmStsSetup.h
CbmMCDataArray::Size
Int_t Size(Int_t fileNumber, Int_t eventNumber)
Definition: CbmMCDataArray.cxx:133
CbmMuchCluster
Data container for MUCH clusters.
Definition: CbmMuchCluster.h:20
CbmMCTrack::GetPdgCode
Int_t GetPdgCode() const
Definition: CbmMCTrack.h:70
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
CbmTrdHit
data class for a reconstructed Energy-4D measurement in the TRD
Definition: CbmTrdHit.h:35
CbmMCDataManager::InitBranch
CbmMCDataArray * InitBranch(const char *name)
Definition: CbmMCDataManager.cxx:106
gMCTracks
static vector< vector< MCTrackDesc > > gMCTracks
Definition: CbmBinnedTrackerQANew.cxx:36
CbmStsSetup::Instance
static CbmStsSetup * Instance()
Definition: CbmStsSetup.cxx:293
CbmGlobalTrack.h
CbmTrdAddress::GetLayerId
static UInt_t GetLayerId(UInt_t address)
Return layer ID from address.
Definition: CbmTrdAddress.h:69
MCTrackDesc::parentInd
Int_t parentInd
Definition: CbmBinnedTrackerQANew.cxx:33
CbmStsHit::GetFrontClusterId
Int_t GetFrontClusterId() const
Definition: CbmStsHit.h:93
CbmHit::GetRefId
Int_t GetRefId() const
Definition: CbmHit.h:72
ECbmModuleId::kTof
@ kTof
Time-of-flight Detector.
CbmStsPoint
Definition: CbmStsPoint.h:27
CbmMatch.h
CbmMuchTrack
Definition: CbmMuchTrack.h:16
CbmStsTrack::GetNofHits
virtual Int_t GetNofHits() const
Definition: CbmStsTrack.h:76
CbmBinnedTrackerQANew
Definition: CbmBinnedTrackerQANew.h:28
CbmCluster::GetNofDigis
Int_t GetNofDigis() const
Number of digis in cluster.
Definition: CbmCluster.h:69
CbmMuchCluster.h
Data container for MUCH clusters.
MCTrackDesc
Definition: CbmBinnedTrackerQANew.cxx:29
CbmTrdCluster
Data Container for TRD clusters.
Definition: CbmTrdCluster.h:23
CbmStsHit
data class for a reconstructed 3-d hit in the STS
Definition: CbmStsHit.h:31
CbmStsDigi.h
CbmGlobalTrack::GetStsTrackIndex
Int_t GetStsTrackIndex() const
Definition: CbmGlobalTrack.h:38
CbmMuchGeoScheme::GetLayerIndex
static Int_t GetLayerIndex(Int_t address)
Definition: CbmMuchGeoScheme.h:71
CbmMuchGeoScheme::GetStationIndex
static Int_t GetStationIndex(Int_t address)
Definition: CbmMuchGeoScheme.h:68
CbmBinnedTrackerQANew.h
CbmMuchPoint.h
tracks
TClonesArray * tracks
Definition: Analyze_matching.h:17
CbmBinnedTrackerQANew::CbmBinnedTrackerQANew
CbmBinnedTrackerQANew()
Definition: CbmBinnedTrackerQANew.cxx:38
CbmStsSetup::GetStationNumber
Int_t GetStationNumber(Int_t address)
Definition: CbmStsSetup.cxx:187
gNofRecoTracks
static int gNofRecoTracks
Definition: CbmBinnedTrackerQANew.cxx:265
CbmTrack::GetHitIndex
Int_t GetHitIndex(Int_t iHit) const
Definition: CbmTrack.h:54
CbmGlobalTrack::GetTrdTrackIndex
Int_t GetTrdTrackIndex() const
Definition: CbmGlobalTrack.h:39
CbmBinnedSettings::Instance
static CbmBinnedSettings * Instance()
Definition: Settings.h:29
MCTrackDesc::mcPointMap
map< pair< ECbmModuleId, Int_t >, pair< UInt_t, set< Int_t > > > mcPointMap
Definition: CbmBinnedTrackerQANew.cxx:31
ClassImp
ClassImp(CbmConverterManager) InitStatus CbmConverterManager
Definition: CbmConverterManager.cxx:12
CbmStsSetup::Init
Bool_t Init(const char *geometryFile=nullptr)
Initialise the setup.
Definition: CbmStsSetup.cxx:201
CbmStsSetup
Class representing the top level of the STS setup.
Definition: CbmStsSetup.h:39
gNofMatchedRecoTracks
static int gNofMatchedRecoTracks
Definition: CbmBinnedTrackerQANew.cxx:266
CbmTrdPoint
Definition: CbmTrdPoint.h:23
CbmMCDataManager
Task class creating and managing CbmMCDataArray objects.
Definition: CbmMCDataManager.h:27
ECbmModuleId::kTrd
@ kTrd
Transition Radiation Detector.
CbmTrdTrack
Definition: CbmTrdTrack.h:22
CbmGlobalTrack
Definition: CbmGlobalTrack.h:26
CbmBinnedTrackerQANew::Init
InitStatus Init()
Definition: CbmBinnedTrackerQANew.cxx:40
CbmStsPoint.h
gNofClones
static int gNofClones
Definition: CbmBinnedTrackerQANew.cxx:267
CbmMCTrack.h
CbmStsHit::GetBackClusterId
Int_t GetBackClusterId() const
Definition: CbmStsHit.h:69
CbmMCTrack
Definition: CbmMCTrack.h:34
CbmTrdPoint::GetModuleAddress
Int_t GetModuleAddress() const
Definition: CbmTrdPoint.h:76
CbmTofPoint.h
m
__m128 m
Definition: L1/vectors/P4_F32vec4.h:26
CbmMuchPixelHit.h
Class for pixel hits in MUCH detector.
CbmBinnedSettings
Definition: Settings.h:27
MCTrackDesc::pdg
Int_t pdg
Definition: CbmBinnedTrackerQANew.cxx:32
ECbmModuleId::kMuch
@ kMuch
Muon detection system.
CbmTofPoint
Geometric intersection of a MC track with a TOFb detector.
Definition: CbmTofPoint.h:40
CbmGlobalTrack::GetTofHitIndex
Int_t GetTofHitIndex() const
Definition: CbmGlobalTrack.h:42
CbmMuchPixelHit
Definition: CbmMuchPixelHit.h:17
CbmTrdCluster.h
Data Container for TRD clusters.
CbmMuchGeoScheme.h
CbmStsTrack::GetStsHitIndex
Int_t GetStsHitIndex(Int_t iHit) const
Definition: CbmStsTrack.h:98
CbmStsCluster.h
Data class for STS clusters.
CbmStsTrack
Definition: CbmStsTrack.h:37
CbmBinnedTrackerQANew::Finish
void Finish()
Definition: CbmBinnedTrackerQANew.cxx:488
ECbmModuleId::kSts
@ kSts
Silicon Tracking System.
CbmCluster::GetDigi
Int_t GetDigi(Int_t index) const
Get digi at position index.
Definition: CbmCluster.h:76
gEventNumber
static Int_t gEventNumber
Definition: CbmBinnedTrackerQANew.cxx:264
CbmStsSetup::IsInit
Bool_t IsInit() const
Initialisation status for sensor parameters.
Definition: CbmStsSetup.h:124
CbmStsHit.h
Data class for a reconstructed hit in the STS.