CbmRoot
CbmTofBuildDigiEvents.cxx
Go to the documentation of this file.
1 
8 
9 
10 #include "CbmMCEventList.h"
11 #include "CbmMatch.h"
12 #include "CbmTimeSlice.h"
13 #include "CbmTofDigi.h"
14 //#include "CbmTofDef.h" TODO
15 
16 #include "FairFileSource.h"
17 #include "FairLogger.h"
18 #include "FairRootManager.h"
19 
20 #include "TClonesArray.h"
21 #include "TMath.h"
22 
23 #include <algorithm>
24 
25 
26 // ---------------------------------------------------------------------------
28  : FairTask("TofBuildDigiEvents")
29  , fFileSource(NULL)
30  , fTimeSliceHeader(NULL)
31  , fTofTimeSliceDigis(NULL)
32  , fDigiMatches(nullptr)
33  , fInputMCEventList(NULL)
34  , fOutputMCEventList(NULL)
35  , fTofEventDigis(NULL)
36  , fdEventWindow(0.)
37  , fNominalTriggerCounterMultiplicity()
38  , fiTriggerMultiplicity(0)
39  , fbPreserveMCBacklinks(kFALSE)
40  , fbMCEventBuilding(kFALSE)
41  , fdEventStartTime(DBL_MIN)
42  , fCounterMultiplicity()
43  , fdIdealEventWindow(1000.)
44  , fProcessedIdealEvents()
45  , fIdealEventStartTimes()
46  , fIdealEventDigis()
47  , fiNEvents(0)
48  , fdDigiToTOffset(0.)
49  , fInactiveCounterSides() {}
50 // ---------------------------------------------------------------------------
51 
52 
53 // ---------------------------------------------------------------------------
55 // ---------------------------------------------------------------------------
56 
57 
58 // ---------------------------------------------------------------------------
59 void CbmTofBuildDigiEvents::Exec(Option_t*) {
60  if (fbMCEventBuilding) {
62 
63  for (Int_t iDigi = 0; iDigi < fTofTimeSliceDigis->GetEntriesFast();
64  iDigi++) {
65  CbmTofDigi* tDigi =
66  dynamic_cast<CbmTofDigi*>(fTofTimeSliceDigis->At(iDigi));
67  CbmMatch* match = dynamic_cast<CbmMatch*>(fDigiMatches->At(iDigi));
68  assert(match);
69 
70  Int_t iDigiAddress = tDigi->GetAddress();
71  Double_t dDigiTime = tDigi->GetTime();
72  Double_t dDigiToT = tDigi->GetTot();
73 
74 
75  Int_t iNMCLinks = match->GetNofLinks();
76 
77  for (Int_t iLink = 0; iLink < iNMCLinks; iLink++) {
78  const CbmLink& tLink = match->GetLink(iLink);
79 
80  // Only consider digis that contain at least one contribution from a
81  // primary source signal. If the digi contains primary signal contributions
82  // from more than one MC event, assign the digi to all MC events found.
83  // TODO: Replace '0' by 'tof::signal_SourcePrimary' upon inclusion of
84  // 'tof/TofTools/CbmTofDef.h' into trunk!
85  if (0 == tLink.GetUniqueID()) {
86  std::pair<Int_t, Int_t> EventID(tLink.GetFile(), tLink.GetEntry());
87 
88  // The MC event is already known.
89  if (fIdealEventStartTimes.find(EventID)
90  != fIdealEventStartTimes.end()) {
91  auto& DigiVector = fIdealEventDigis.at(EventID);
92 
94  // deep copy construction including 'CbmDigi::fMatch'
95  DigiVector.push_back(new CbmTofDigi(*tDigi));
96  } else {
97  // shallow construction excluding 'CbmDigi::fMatch'
98  DigiVector.push_back(
99  new CbmTofDigi(iDigiAddress, dDigiTime, dDigiToT));
100  }
101  }
102  // The MC event is not known yet.
103  else {
104  // Make sure that a late digi from an event that has already been
105  // processed and written to disk (i.e. the time difference to the
106  // earliest digi in the same event is larger than 'fdIdealEventWindow')
107  // does not trigger separate event processing for itself only
108  // (and possibly a few additional latecomers).
109  if (fProcessedIdealEvents.find(EventID)
110  == fProcessedIdealEvents.end()) {
111  fIdealEventStartTimes.emplace(EventID, dDigiTime);
112  fIdealEventDigis.emplace(EventID, std::vector<CbmTofDigi*>());
113 
114  auto& DigiVector = fIdealEventDigis.at(EventID);
115 
116  if (fbPreserveMCBacklinks) {
117  // deep copy construction including 'CbmDigi::fMatch'
118  DigiVector.push_back(new CbmTofDigi(*tDigi));
119  } else {
120  // shallow construction excluding 'CbmDigi::fMatch'
121  DigiVector.push_back(
122  new CbmTofDigi(iDigiAddress, dDigiTime, dDigiToT));
123  }
124  }
125  }
126  }
127  }
128  }
129  } else {
130  for (Int_t iDigi = 0; iDigi < fTofTimeSliceDigis->GetEntriesFast();
131  iDigi++) {
132  CbmTofDigi* tDigi =
133  dynamic_cast<CbmTofDigi*>(fTofTimeSliceDigis->At(iDigi));
134 
135  Int_t iDigiModuleType = tDigi->GetType();
136  Int_t iDigiModuleIndex = tDigi->GetSm();
137  Int_t iDigiCounterIndex = tDigi->GetRpc();
138  Int_t iDigiCounterSide = tDigi->GetSide();
139  Int_t iDigiAddress = tDigi->GetAddress();
140  Double_t dDigiTime = tDigi->GetTime();
141  Double_t dDigiToT = tDigi->GetTot();
142 
143 
144  if (dDigiTime - fdEventStartTime > fdEventWindow) {
145  std::map<std::tuple<Int_t, Int_t, Int_t>, UChar_t>
146  ActualTriggerCounterMultiplicity;
147 
148  std::set_intersection(
149  fCounterMultiplicity.begin(),
150  fCounterMultiplicity.end(),
153  std::inserter(ActualTriggerCounterMultiplicity,
154  ActualTriggerCounterMultiplicity.begin()));
155 
156  if (ActualTriggerCounterMultiplicity.size()
157  >= static_cast<size_t>(fiTriggerMultiplicity)) {
159 
160  FairRootManager::Instance()->Fill();
161  fiNEvents++;
163  fTofEventDigis->Delete();
164  } else {
165  fTofEventDigis->Delete();
166  }
167 
168  fCounterMultiplicity.clear();
169 
170  fdEventStartTime = dDigiTime;
171  }
172 
173 
174  fCounterMultiplicity[std::make_tuple(
175  iDigiModuleType, iDigiModuleIndex, iDigiCounterIndex)] |=
176  1 << iDigiCounterSide;
177 
178  auto CounterSideTuple = std::make_tuple(
179  iDigiModuleType, iDigiModuleIndex, iDigiCounterIndex, iDigiCounterSide);
180 
181  if (fInactiveCounterSides.find(CounterSideTuple)
182  == fInactiveCounterSides.end()) {
183  CbmTofDigi* tEventDigi(NULL);
184 
185  if (fbPreserveMCBacklinks) {
186  // deep copy construction including 'CbmDigi::fMatch'
187  tEventDigi = new ((*fTofEventDigis)[fTofEventDigis->GetEntriesFast()])
188  CbmTofDigi(*tDigi);
189  } else {
190  // shallow construction excluding 'CbmDigi::fMatch'
191  tEventDigi = new ((*fTofEventDigis)[fTofEventDigis->GetEntriesFast()])
192  CbmTofDigi(iDigiAddress, dDigiTime, dDigiToT);
193  }
194 
195  tEventDigi->SetTot(tEventDigi->GetTot() + fdDigiToTOffset);
196  }
197  }
198  }
199 }
200 // ---------------------------------------------------------------------------
201 
202 
203 // ---------------------------------------------------------------------------
205  if (!FairRootManager::Instance()) {
206  LOG(error) << "FairRootManager not found.";
207  return kFATAL;
208  }
209 
210  fFileSource =
211  dynamic_cast<FairFileSource*>(FairRootManager::Instance()->GetSource());
212  if (!fFileSource) {
213  LOG(error) << "Could not get pointer to FairFileSource.";
214  return kFATAL;
215  }
216 
217  fTimeSliceHeader = dynamic_cast<CbmTimeSlice*>(
218  FairRootManager::Instance()->GetObject("TimeSlice."));
219  if (!fTimeSliceHeader) {
220  LOG(error)
221  << "Could not retrieve branch 'TimeSlice.' from FairRootManager.";
222  return kFATAL;
223  }
224 
225  fTofTimeSliceDigis = dynamic_cast<TClonesArray*>(
226  FairRootManager::Instance()->GetObject("TofDigiExp"));
227  if (!fTofTimeSliceDigis) {
228  LOG(error)
229  << "Could not retrieve branch 'TofDigiExp' from FairRootManager.";
230  return kFATAL;
231  }
232 
233  fDigiMatches = dynamic_cast<TClonesArray*>(
234  FairRootManager::Instance()->GetObject("TofDigiMatch"));
235  if (!fDigiMatches) {
236  LOG(error)
237  << "Could not retrieve branch 'TofDigiMatch' from FairRootManager.";
238  return kFATAL;
239  }
240 
241  fInputMCEventList = dynamic_cast<CbmMCEventList*>(
242  FairRootManager::Instance()->GetObject("MCEventList."));
243  if (!fInputMCEventList) {
244  LOG(error)
245  << "Could not retrieve branch 'MCEventList.' from FairRootManager.";
246  return kFATAL;
247  }
248 
249  if (FairRootManager::Instance()->GetObject("TofPointTB")) {
250  LOG(error) << "Timeslice branch with MC points found. Event building would "
251  "not work properly.";
252  return kFATAL;
253  }
254 
255 
257  FairRootManager::Instance()->Register("EventList.",
258  "EventList",
260  IsOutputBranchPersistent("EventList."));
261 
262  fTofEventDigis = new TClonesArray("CbmTofDigi", 100);
263  FairRootManager::Instance()->Register("CbmTofDigi",
264  "TOF event digis",
266  IsOutputBranchPersistent("CbmTofDigi"));
267 
268 
269  if (0. >= fdEventWindow) { fbMCEventBuilding = kTRUE; }
270 
272 
274  < static_cast<size_t>(fiTriggerMultiplicity)) {
276  }
277 
278  return kSUCCESS;
279 }
280 // ---------------------------------------------------------------------------
281 
282 
283 // ---------------------------------------------------------------------------
285  if (fbMCEventBuilding) {
286  // With O(s) of off-spill noise (not eligible for MC event building) stored
287  // in several timeslices (the processing of each causing an 'Exec' call by
288  // the framework) following the final spill there should not be any digis
289  // related to MC events left for processing at this point.
290  ProcessIdealEvents(DBL_MAX);
291  } else {
292  // The remaining digis in the buffer do not cover a time interval of
293  // 'fdEventWindow' and, in consequence, do not qualify for event building.
294  fTofEventDigis->Delete();
295  fCounterMultiplicity.clear();
296  }
297 }
298 // ---------------------------------------------------------------------------
299 
300 
301 // ---------------------------------------------------------------------------
303  Int_t iModuleIndex,
304  Int_t iCounterIndex,
305  Int_t iNCounterSides) {
307  std::make_tuple(iModuleType, iModuleIndex, iCounterIndex),
308  (1 == iNCounterSides) ? 1 : 3);
309 }
310 // ---------------------------------------------------------------------------
311 
312 
313 // ---------------------------------------------------------------------------
314 void CbmTofBuildDigiEvents::ProcessIdealEvents(Double_t dProcessingTime) {
315  for (auto itEvent = fIdealEventStartTimes.cbegin();
316  itEvent != fIdealEventStartTimes.cend();) {
317  auto EventID = itEvent->first;
318  Double_t dEventStartTime = itEvent->second;
319 
320  if (dProcessingTime - dEventStartTime > fdIdealEventWindow) {
321  for (auto& tDigi : fIdealEventDigis.at(EventID)) {
322  Int_t iDigiModuleType = tDigi->GetType();
323  Int_t iDigiModuleIndex = tDigi->GetSm();
324  Int_t iDigiCounterIndex = tDigi->GetRpc();
325  Int_t iDigiCounterSide = tDigi->GetSide();
326 
327  fCounterMultiplicity[std::make_tuple(
328  iDigiModuleType, iDigiModuleIndex, iDigiCounterIndex)] |=
329  1 << iDigiCounterSide;
330 
331  auto CounterSideTuple = std::make_tuple(iDigiModuleType,
332  iDigiModuleIndex,
333  iDigiCounterIndex,
334  iDigiCounterSide);
335 
336  if (fInactiveCounterSides.find(CounterSideTuple)
337  == fInactiveCounterSides.end()) {
338  // deep copy construction including 'CbmDigi::fMatch' (only if already deep-copied in 'Exec')
339  CbmTofDigi* tEventDigi =
340  new ((*fTofEventDigis)[fTofEventDigis->GetEntriesFast()])
341  CbmTofDigi(*tDigi);
342  tEventDigi->SetTot(tEventDigi->GetTot() + fdDigiToTOffset);
343  }
344 
345  delete tDigi;
346  }
347 
348  std::map<std::tuple<Int_t, Int_t, Int_t>, UChar_t>
349  ActualTriggerCounterMultiplicity;
350 
351  std::set_intersection(
352  fCounterMultiplicity.begin(),
353  fCounterMultiplicity.end(),
356  std::inserter(ActualTriggerCounterMultiplicity,
357  ActualTriggerCounterMultiplicity.begin()));
358 
359  if (ActualTriggerCounterMultiplicity.size()
360  >= static_cast<size_t>(fiTriggerMultiplicity)) {
362 
363  FairRootManager::Instance()->Fill();
364  fiNEvents++;
366  fTofEventDigis->Delete();
367  } else {
368  fTofEventDigis->Delete();
369  }
370 
371  fCounterMultiplicity.clear();
372  fIdealEventDigis.erase(EventID);
373  fProcessedIdealEvents.emplace(EventID);
374 
375  itEvent = fIdealEventStartTimes.erase(itEvent);
376  } else {
377  ++itEvent;
378  }
379  }
380 }
381 // ---------------------------------------------------------------------------
382 
383 
384 // ---------------------------------------------------------------------------
386  std::set<std::pair<Int_t, Int_t>> MCEventSet;
387 
388  for (Int_t iDigi = 0; iDigi < fTofEventDigis->GetEntriesFast(); iDigi++) {
389  //CbmTofDigi* tDigi = dynamic_cast<CbmTofDigi*>(fTofEventDigis->At(iDigi)); (VF) not used
390  CbmMatch* match = dynamic_cast<CbmMatch*>(fDigiMatches->At(iDigi));
391  Int_t iNMCLinks = match->GetNofLinks();
392 
393  for (Int_t iLink = 0; iLink < iNMCLinks; iLink++) {
394  const CbmLink& tLink = match->GetLink(iLink);
395 
396  Int_t iFileIndex = tLink.GetFile();
397  Int_t iEventIndex = tLink.GetEntry();
398 
399  // Collect original MC event affiliations of digis attributed to
400  // the current reconstructed event.
401  if (-1 < iFileIndex && -1 < iEventIndex) {
402  MCEventSet.emplace(iFileIndex, iEventIndex);
403  }
404  }
405  }
406 
407  // Read the respective start times of the original MC events contributing to
408  // the reconstructed event from the input MC event list and make them
409  // available in the output MC event list.
410  for (auto const& MCEvent : MCEventSet) {
411  Int_t iFileIndex = MCEvent.first;
412  Int_t iEventIndex = MCEvent.second;
413 
414  Double_t dStartTime =
415  fInputMCEventList->GetEventTime(iEventIndex, iFileIndex);
416 
417  if (-1. != dStartTime) {
418  fOutputMCEventList->Insert(iEventIndex, iFileIndex, dStartTime);
419  } else {
420  LOG(fatal) << Form(
421  "Could not find MC event (%d, %d) in the input MC event list.",
422  iFileIndex,
423  iEventIndex);
424  }
425  }
426 
428 }
429 // ---------------------------------------------------------------------------
430 
431 
432 // ---------------------------------------------------------------------------
434  Int_t iModuleIndex,
435  Int_t iCounterIndex,
436  Int_t iCounterSide) {
437  fInactiveCounterSides.emplace(
438  std::make_tuple(iModuleType, iModuleIndex, iCounterIndex, iCounterSide));
439 }
440 // ---------------------------------------------------------------------------
441 
442 
CbmTofBuildDigiEvents
...
Definition: CbmTofBuildDigiEvents.h:32
CbmMatch
Definition: CbmMatch.h:22
CbmMCEventList::Sort
void Sort()
Sort the list.
Definition: CbmMCEventList.cxx:133
CbmTofBuildDigiEvents::fNominalTriggerCounterMultiplicity
std::map< std::tuple< Int_t, Int_t, Int_t >, UChar_t > fNominalTriggerCounterMultiplicity
Definition: CbmTofBuildDigiEvents.h:85
CbmTofBuildDigiEvents::SetIgnoreCounterSide
void SetIgnoreCounterSide(Int_t iModuleType, Int_t iModuleIndex, Int_t iCounterIndex, Int_t iCounterSide)
Definition: CbmTofBuildDigiEvents.cxx:433
CbmMCEventList::Insert
Bool_t Insert(UInt_t event, UInt_t file, Double_t time)
Definition: CbmMCEventList.cxx:116
CbmMatch::GetLink
const CbmLink & GetLink(Int_t i) const
Definition: CbmMatch.h:35
CbmMatch::GetNofLinks
Int_t GetNofLinks() const
Definition: CbmMatch.h:38
CbmTofBuildDigiEvents::fDigiMatches
TClonesArray * fDigiMatches
Definition: CbmTofBuildDigiEvents.h:79
CbmTofBuildDigiEvents::fProcessedIdealEvents
std::set< std::pair< Int_t, Int_t > > fProcessedIdealEvents
Definition: CbmTofBuildDigiEvents.h:92
CbmTofBuildDigiEvents::Exec
virtual void Exec(Option_t *option)
Definition: CbmTofBuildDigiEvents.cxx:59
CbmTofBuildDigiEvents::FillMCEventList
void FillMCEventList()
Definition: CbmTofBuildDigiEvents.cxx:385
CbmTofBuildDigiEvents.h
CbmTofDigi::SetTot
void SetTot(Double_t tot)
Definition: CbmTofDigi.h:154
CbmTofDigi.h
CbmTimeSlice::GetStartTime
Double_t GetStartTime() const
Definition: CbmTimeSlice.h:108
CbmMatch.h
CbmTofDigi::GetSm
Double_t GetSm() const
Sm.
Definition: CbmTofDigi.h:124
CbmTimeSlice.h
CbmTofBuildDigiEvents::Init
virtual InitStatus Init()
Definition: CbmTofBuildDigiEvents.cxx:204
CbmTofBuildDigiEvents::~CbmTofBuildDigiEvents
virtual ~CbmTofBuildDigiEvents()
Definition: CbmTofBuildDigiEvents.cxx:54
CbmTofBuildDigiEvents::SetTriggerCounter
void SetTriggerCounter(Int_t iModuleType, Int_t iModuleIndex, Int_t iCounterIndex, Int_t iNCounterSides)
Definition: CbmTofBuildDigiEvents.cxx:302
CbmTofDigi::GetSide
Double_t GetSide() const
Channel Side.
Definition: CbmTofDigi.h:142
CbmTofBuildDigiEvents::fInputMCEventList
CbmMCEventList * fInputMCEventList
Definition: CbmTofBuildDigiEvents.h:80
CbmTofDigi::GetType
Double_t GetType() const
Sm Type .
Definition: CbmTofDigi.h:128
CbmTofBuildDigiEvents::fTofEventDigis
TClonesArray * fTofEventDigis
Definition: CbmTofBuildDigiEvents.h:82
CbmTofBuildDigiEvents::fIdealEventStartTimes
std::map< std::pair< Int_t, Int_t >, Double_t > fIdealEventStartTimes
Definition: CbmTofBuildDigiEvents.h:93
CbmMCEventList::GetEventTime
Double_t GetEventTime(UInt_t event, UInt_t file)
Event start time.
Definition: CbmMCEventList.cxx:86
CbmTofBuildDigiEvents::fdEventStartTime
Double_t fdEventStartTime
Definition: CbmTofBuildDigiEvents.h:89
CbmTofBuildDigiEvents::fdIdealEventWindow
Double_t fdIdealEventWindow
Definition: CbmTofBuildDigiEvents.h:91
CbmTofDigi::GetTime
Double_t GetTime() const
Inherited from CbmDigi.
Definition: CbmTofDigi.h:111
CbmTofBuildDigiEvents::fFileSource
FairFileSource * fFileSource
Definition: CbmTofBuildDigiEvents.h:76
ClassImp
ClassImp(CbmConverterManager) InitStatus CbmConverterManager
Definition: CbmConverterManager.cxx:12
CbmMCEventList
Container class for MC events with number, file and start time.
Definition: CbmMCEventList.h:38
CbmTofBuildDigiEvents::ProcessIdealEvents
void ProcessIdealEvents(Double_t dProcessingTime)
Definition: CbmTofBuildDigiEvents.cxx:314
CbmTofBuildDigiEvents::fTofTimeSliceDigis
TClonesArray * fTofTimeSliceDigis
Definition: CbmTofBuildDigiEvents.h:78
CbmTofBuildDigiEvents::CbmTofBuildDigiEvents
CbmTofBuildDigiEvents()
Definition: CbmTofBuildDigiEvents.cxx:27
CbmTofDigi
Data class for expanded digital TOF information.
Definition: CbmTofDigi.h:38
CbmTofBuildDigiEvents::fOutputMCEventList
CbmMCEventList * fOutputMCEventList
Definition: CbmTofBuildDigiEvents.h:81
CbmMCEventList.h
CbmTofDigi::GetRpc
Double_t GetRpc() const
Detector aka Module aka RPC .
Definition: CbmTofDigi.h:132
CbmTofBuildDigiEvents::fTimeSliceHeader
CbmTimeSlice * fTimeSliceHeader
Definition: CbmTofBuildDigiEvents.h:77
CbmTofBuildDigiEvents::Finish
virtual void Finish()
Definition: CbmTofBuildDigiEvents.cxx:284
CbmTofBuildDigiEvents::fiNEvents
Int_t fiNEvents
Definition: CbmTofBuildDigiEvents.h:95
CbmTimeSlice
Bookkeeping of time-slice content.
Definition: CbmTimeSlice.h:29
CbmTofBuildDigiEvents::fdDigiToTOffset
Double_t fdDigiToTOffset
Definition: CbmTofBuildDigiEvents.h:96
CbmTofBuildDigiEvents::fCounterMultiplicity
std::map< std::tuple< Int_t, Int_t, Int_t >, UChar_t > fCounterMultiplicity
Definition: CbmTofBuildDigiEvents.h:90
CbmMCEventList::Clear
virtual void Clear(Option_t *)
Delete all event entries.
Definition: CbmMCEventList.h:50
CbmTofBuildDigiEvents::fIdealEventDigis
std::map< std::pair< Int_t, Int_t >, std::vector< CbmTofDigi * > > fIdealEventDigis
Definition: CbmTofBuildDigiEvents.h:94
CbmTofBuildDigiEvents::fbMCEventBuilding
Bool_t fbMCEventBuilding
Definition: CbmTofBuildDigiEvents.h:88
CbmTofBuildDigiEvents::fbPreserveMCBacklinks
Bool_t fbPreserveMCBacklinks
Definition: CbmTofBuildDigiEvents.h:87
CbmTofDigi::GetTot
Double_t GetTot() const
Alias for GetCharge.
Definition: CbmTofDigi.h:120
CbmTofBuildDigiEvents::fiTriggerMultiplicity
Int_t fiTriggerMultiplicity
Definition: CbmTofBuildDigiEvents.h:86
CbmTofDigi::GetAddress
Int_t GetAddress() const
Inherited from CbmDigi.
Definition: CbmTofDigi.h:98
CbmTofBuildDigiEvents::fInactiveCounterSides
std::set< std::tuple< Int_t, Int_t, Int_t, Int_t > > fInactiveCounterSides
Definition: CbmTofBuildDigiEvents.h:97
CbmTofBuildDigiEvents::fdEventWindow
Double_t fdEventWindow
Definition: CbmTofBuildDigiEvents.h:83