ATLAS Offline Software
ALFA_PileUpTool.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 #include "ALFA_Digitization/ALFA_PileUpTool.h"
6 
7 #include "ALFA_SimEv/ALFA_HitCollection.h"
8 #include "ALFA_SimEv/ALFA_Hit.h"
9 #include "ALFA_SimEv/ALFA_ODHitCollection.h"
10 
11 #include "ALFA_RawEv/ALFA_DigitCollection.h"
12 #include "ALFA_SimEv/ALFA_ODHit.h"
13 #include "ALFA_RawEv/ALFA_ODDigitCollection.h"
14 
15 #include "AthenaKernel/RNGWrapper.h"
16 #include "CLHEP/Random/RandomEngine.h"
17 #include "CLHEP/Random/RandFlat.h"
18 #include "CLHEP/Random/RandGaussZiggurat.h"
19 #include "CLHEP/Random/RandPoissonQ.h"
20 
21 #include <algorithm>
22 
23 #include <map>
24 
25 #include "Identifier/Identifier.h"
26 #include "AthenaBaseComps/AthMsgStreamMacros.h"
27 
28 // from now old part
29 
30 #include <cmath>
31 #include <functional>
32 
33 //for truth particles
34 #include "AtlasHepMC/GenParticle.h"
35 #include "GeneratorObjects/McEventCollection.h"
36 #include "AtlasHepMC/GenEvent.h"
37 #include "AtlasHepMC/GenVertex.h"
38 #include "CLHEP/Vector/LorentzVector.h"
39 #include "AtlasHepMC/SimpleVector.h"
40 
41 #include "TTree.h"
42 
43 #include "GaudiKernel/ITHistSvc.h"
44 
45 #include "PathResolver/PathResolver.h"
46 
47 // *********************************************************************************
48 
49 ALFA_PileUpTool::ALFA_PileUpTool(const std::string& type,
50  const std::string& name,
51  const IInterface* parent) :
52  PileUpToolBase(type, name, parent),
53  m_SimHitCollectionName ("ALFA_HitCollection"),
54  m_SimODHitCollectionName ("ALFA_ODHitCollection"),
55  m_key_DigitCollection ("ALFA_DigitCollection"),
56  m_key_ODDigitCollection ("ALFA_ODDigitCollection"),
57  m_digitCollection (nullptr),
58  m_ODdigitCollection (nullptr),
59  m_sigma0 (0.1),
60  m_sigma1 (1.0),
61  m_meanE_dep (0.0863), //MeV
62  m_meanN_photo (4.11),
63  m_mean (0.0),
64  m_stdDev (1.0),
65  m_AmplitudeCut (0.24),
66  m_mergedALFA_HitList (nullptr),
67  m_mergedALFA_ODHitList (nullptr)
68 {
69  // Properties
70  declareProperty("SimHitCollection", m_SimHitCollectionName, "Name of the input ALFA Sim Hit Collection of simulated hits");
71  declareProperty("SimODHitCollection", m_SimODHitCollectionName, "Name of the input ALFA Sim OD Hit Collection of simulated hits");
72 
73  declareProperty("ALFA_DigitCollection", m_key_DigitCollection, "Name of the Collection to hold the output from the ALFA main detector digitization");
74  declareProperty("ALFA_ODDigitCollection", m_key_ODDigitCollection, "Name of the Collection to hold the output from the ALFA OD digitization");
75 
76 
77 // declareProperty("fillRootTree", m_fillRootTree);
78 
79 
80  declareProperty("sigma0", m_sigma0);
81  declareProperty("sigma1", m_sigma1);
82  declareProperty("meanE_dep", m_meanE_dep);
83  declareProperty("meanN_photo", m_meanN_photo);
84  declareProperty("mean", m_mean);
85  declareProperty("stdDev", m_stdDev);
86  declareProperty("AmplitudeCut", m_AmplitudeCut);
87 // declareProperty("fillRootTree", m_fillRootTree);
88 }
89 
90 
91 
92 
93 StatusCode ALFA_PileUpTool::initialize(){
94 
95  ATH_MSG_DEBUG ("ALFA_digiti initialize() called");
96 
97  ATH_MSG_DEBUG (" Pedestal: " << m_AmplitudeCut);
98 
99  ATH_CHECK (m_randomSvc.retrieve());
100  ATH_MSG_DEBUG ( "Retrieved RandomNumber Service" );
101 
102  ATH_CHECK (m_mergeSvc.retrieve());
103  ATH_MSG_DEBUG("Retrieved PileUpMergeSvc");
104 
105  m_mergedALFA_HitList = new ALFA_HitCollection();
106  m_mergedALFA_ODHitList = new ALFA_ODHitCollection();
107 
108  return StatusCode::SUCCESS;
109 }
110 
111 
112 
113 
114 // *********************************************************************************
115 
116 StatusCode ALFA_PileUpTool::processAllSubEvents(const EventContext& ctx) {
117 
118  ATH_MSG_DEBUG ("ALFA_PileUpTool::processAllSubEvents()");
119 
120  using TimedALFAHitCollList = PileUpMergeSvc::TimedList<ALFA_HitCollection>::type;
121  using TimedALFAODHitCollList = PileUpMergeSvc::TimedList<ALFA_ODHitCollection>::type;
122 
123  TimedALFAHitCollList tHitCollList;
124  TimedALFAODHitCollList tODHitCollList;
125 
126 
127  if (m_mergeSvc->retrieveSubEvtsData(m_SimHitCollectionName, tHitCollList).isFailure()) {
128 
129  ATH_MSG_FATAL ( "Could not fill TimedALFAHitCollList" ); return StatusCode::FAILURE;
130  }
131  else { ATH_MSG_DEBUG ( "Retrieved TimedALFAHitCollList" ); }
132 
133  ATH_MSG_DEBUG ( "PileUp: Merge " << tHitCollList.size() << " ALFA_HitCollection with key" << m_SimHitCollectionName );
134 
135 
136  if (m_mergeSvc->retrieveSubEvtsData(m_SimODHitCollectionName, tODHitCollList).isFailure()) {
137 
138  ATH_MSG_FATAL ( "Could not fill TimedALFAODHitCollList" ); return StatusCode::FAILURE;
139  }
140  else { ATH_MSG_DEBUG ( "Retrieved TimedALFAODHitCollList" ); }
141 
142  ATH_MSG_DEBUG ( "PileUp: Merge " << tODHitCollList.size() << " ALFA_ODHitCollection with key " << m_SimODHitCollectionName );
143 
144 
145  TimedHitCollection<ALFA_Hit> tALFAhit;
146 
147  TimedALFAHitCollList::iterator iHitColl (tHitCollList.begin());
148  TimedALFAHitCollList::iterator eHitColl (tHitCollList.end());
149 
150  while (iHitColl != eHitColl) {
151 
152  const ALFA_HitCollection* tmpColl(iHitColl->second);
153 
154  tALFAhit.insert(iHitColl->first, tmpColl);
155 
156  ATH_MSG_DEBUG ( " ALFA_HitCollection found with " << tmpColl->size() << " hits " << iHitColl->first );
157 
158  ++iHitColl;
159  }
160 
161  TimedHitCollection<ALFA_ODHit> tALFAODhit;
162 
163  TimedALFAODHitCollList::iterator iODHitColl (tODHitCollList.begin());
164  TimedALFAODHitCollList::iterator eODHitColl (tODHitCollList.end());
165 
166  while (iODHitColl != eODHitColl) {
167 
168  const ALFA_ODHitCollection* tmpColl(iODHitColl->second);
169 
170  tALFAODhit.insert(iODHitColl->first, tmpColl);
171 
172  ATH_MSG_DEBUG ( " ALFA_ODHitCollection found with " << tmpColl->size() << " hits " << iODHitColl->first );
173 
174  ++iODHitColl;
175  }
176 
177  // preparing containers
178 
179  if (recordCollection(this->evtStore(), m_key_DigitCollection).isFailure()) {
180 
181  ATH_MSG_FATAL ( " ALFA_PileUpTool::processAllSubEvents(): Could not record the empty ALFA digit container in StoreGate " ); return StatusCode::FAILURE;
182  }
183  else { ATH_MSG_DEBUG ( " ALFA_PileUpTool::processAllSubEvents(): ALFA Digit container is recorded in StoreGate " ); }
184 
185 
186  if (recordODCollection(this->evtStore(), m_key_ODDigitCollection).isFailure()) {
187 
188  ATH_MSG_FATAL ( " ALFA_PileUpTool::processAllSubEvents(): Could not record the empty ALFA OD digit container in StoreGate " ); return StatusCode::FAILURE;
189  }
190  else { ATH_MSG_DEBUG ( " ALFA_PileUpTool::processAllSubEvents(): ALFA OD Digit container is recorded in StoreGate " ); }
191 
192 
193  // filling containers
194 
195  ALFA_MD_info(tALFAhit);
196 
197  // Prepare RNG Service
198  ATHRNG::RNGWrapper* rngWrapper = m_randomSvc->getEngine(this, m_randomStreamName);
199  rngWrapper->setSeed( m_randomStreamName, ctx );
200  CLHEP::HepRandomEngine* rngEngine = rngWrapper->getEngine(ctx);
201  StatusCode sc = fill_MD_DigitCollection(rngEngine);
202  if (sc.isFailure()) {
203  ATH_MSG_WARNING ("ALFA_PileUpTool::fill_MD_DigitCollection() failed");
204  return StatusCode::SUCCESS;
205  }
206 
207 
208  ALFA_OD_info(tALFAODhit);
209 
210  sc = fill_OD_DigitCollection(rngEngine);
211  if (sc.isFailure()) {
212  ATH_MSG_WARNING ("ALFA_PileUpTool::fill_OD_DigitCollection() failed");
213  return StatusCode::SUCCESS;
214  }
215 
216  return StatusCode::SUCCESS;
217 }
218 
219 
220 
221 
222 StatusCode ALFA_PileUpTool::prepareEvent(const EventContext& /*ctx*/, const unsigned int nInputEvents){
223 
224  ATH_MSG_DEBUG ( "ALFA_PileUpTool::prepareEvent() called for " << nInputEvents << " input events" );
225 
226  StatusCode sc = recordCollection(this->evtStore(), m_key_DigitCollection);
227 
228  if (sc.isFailure()) { ATH_MSG_FATAL ( " ALFA_PileUpTool::prepareEvent(): Could not record the empty digit container in StoreGate " ); return sc; }
229  else { ATH_MSG_DEBUG ( " ALFA_PileUpTool::prepareEvent(): Digit container is recorded in StoreGate " ); }
230 
231  sc = recordODCollection(this->evtStore(), m_key_ODDigitCollection);
232 
233  if (sc.isFailure()) { ATH_MSG_FATAL ( " ALFA_PileUpTool::prepareEvent(): Could not record the empty digit OD container in StoreGate " ); return sc; }
234  else { ATH_MSG_DEBUG ( " ALFA_PileUpTool::prepareEvent(): Digit OD container is recorded in StoreGate " ); }
235 
236  m_mergedALFA_HitList->clear();
237  m_mergedALFA_ODHitList->clear();
238 
239  return StatusCode::SUCCESS;
240 }
241 
242 
243 StatusCode ALFA_PileUpTool::processBunchXing(int bunchXing,
244  SubEventIterator bSubEvents,
245  SubEventIterator eSubEvents) {
246  ATH_MSG_DEBUG ( "ALFA_PileUpTool::processBunchXing() " << bunchXing );
247  SubEventIterator iEvt = bSubEvents;
248  for (; iEvt!=eSubEvents; ++iEvt) {
249  StoreGateSvc& seStore = *iEvt->ptr()->evtStore();
250  //PileUpTimeEventIndex thisEventIndex = PileUpTimeEventIndex(static_cast<int>(iEvt->time()),iEvt->index()); not used
251  ATH_MSG_VERBOSE("SubEvt StoreGate " << seStore.name() << " :"
252  << " bunch crossing : " << bunchXing
253  << " time offset : " << iEvt->time()
254  << " event number : " << iEvt->ptr()->eventNumber()
255  << " run number : " << iEvt->ptr()->runNumber()
256  );
257 
258 
259  const ALFA_HitCollection* tmpHitColl = nullptr;
260 
261  if (!seStore.retrieve(tmpHitColl, m_SimHitCollectionName).isSuccess()) {
262 
263  ATH_MSG_ERROR ("SubEvent ALFA_HitCollection not found in StoreGate " << seStore.name() );
264 
265  return StatusCode::FAILURE;
266  }
267 
268  ATH_MSG_DEBUG ("SubEvent, ALFA_HitCollection found with " << tmpHitColl->size() << " hits");
269 
270  ALFA_HitCollection::const_iterator iHitColl = tmpHitColl->begin();
271  ALFA_HitCollection::const_iterator eHitColl = tmpHitColl->end();
272 
273  for (; iHitColl!=eHitColl; ++iHitColl) m_mergedALFA_HitList->push_back((*iHitColl));
274 
275  const ALFA_ODHitCollection* tmpODHitColl = nullptr;
276 
277  if (!seStore.retrieve(tmpODHitColl, m_SimODHitCollectionName).isSuccess()) {
278 
279  ATH_MSG_ERROR ( "SubEvent ALFA_ODHitCollection not found in StoreGate " << seStore.name() );
280 
281  return StatusCode::FAILURE;
282  }
283 
284  ATH_MSG_DEBUG ( "ALFA_ODHitCollection found with " << tmpODHitColl->size() << " hits" );
285 
286  ALFA_ODHitCollection::const_iterator iODHitColl = tmpODHitColl->begin();
287  ALFA_ODHitCollection::const_iterator eODHitColl = tmpODHitColl->end();
288 
289  for (; iODHitColl!=eODHitColl; ++iODHitColl) m_mergedALFA_ODHitList->push_back((*iODHitColl));
290  }
291 
292  return StatusCode::SUCCESS;
293 }
294 
295 
296 
297 
298 StatusCode ALFA_PileUpTool::mergeEvent(const EventContext& ctx){
299 
300 
301  ALFA_MD_info(m_mergedALFA_HitList);
302 
303  // Prepare RNG Service
304  ATHRNG::RNGWrapper* rngWrapper = m_randomSvc->getEngine(this, m_randomStreamName);
305  rngWrapper->setSeed( m_randomStreamName, ctx );
306  CLHEP::HepRandomEngine* rngEngine = rngWrapper->getEngine(ctx);
307  StatusCode sc = fill_MD_DigitCollection(rngEngine);
308  if (sc.isFailure()) {
309  ATH_MSG_WARNING ("ALFA_PileUpTool::fill_MD_DigitCollection() failed");
310  return StatusCode::SUCCESS;
311  }
312 
313 
314  ALFA_OD_info(m_mergedALFA_ODHitList);
315 
316  sc = fill_OD_DigitCollection(rngEngine);
317  if (sc.isFailure()) {
318  ATH_MSG_WARNING ("ALFA_PileUpTool::fill_OD_DigitCollection() failed");
319  return StatusCode::SUCCESS;
320  }
321 
322  return StatusCode::SUCCESS;
323 }
324 
325 StatusCode ALFA_PileUpTool::finalize() { return StatusCode::SUCCESS; }
326 
327 
328 
329 StatusCode ALFA_PileUpTool::recordCollection(ServiceHandle<StoreGateSvc>& evtStore, const std::string& key_digitCnt)
330 {
331  ATH_MSG_DEBUG ("ALFA_Digitization::recordCollection");
332 
333  m_digitCollection = new ALFA_DigitCollection();
334 
335  StatusCode sc = evtStore->record(m_digitCollection, key_digitCnt);
336 
337  return sc;
338 }
339 
340 
341 StatusCode ALFA_PileUpTool::recordODCollection(ServiceHandle<StoreGateSvc>& evtStore, const std::string& key_ODdigitCnt) {
342 
343  ATH_MSG_DEBUG ("ALFA_Digitization::recordODCollection");
344 
345  m_ODdigitCollection = new ALFA_ODDigitCollection();
346 
347  StatusCode sc = evtStore->record(m_ODdigitCollection, key_ODdigitCnt);
348 
349  return sc;
350 }
351 
352 
353 
354 
355 void ALFA_PileUpTool::ALFA_MD_info(TimedHitCollection<ALFA_Hit>& tHitALFA)
356 {
357 
358 // cleaning
359 
360  for ( int l = 0; l < 8; l++ )
361  {
362  for ( int i = 0; i < 20; i++ )
363  {
364  for ( int j = 0; j < 64; j++ )
365  {
366 
367  m_E_fib[l][i][j] = 0.;
368  m_fibres[l][i][j] = 0;
369 
370  }
371  }
372  }
373 
374 // cleaning end
375 
376  int fiber, plate, sign, station;
377 
378  TimedHitCollection<ALFA_Hit> thpc = tHitALFA;
379  TimedHitCollection<ALFA_Hit>::const_iterator i, e, it;
380 
381  while (thpc.nextDetectorElement(i, e)) for (it = i; it != e; ++it) {
382 
383  station = (*it)->GetStationNumber();
384  plate = (*it)->GetPlateNumber();
385  fiber = (*it)->GetFiberNumber();
386  sign = (*it)->GetSignFiber();
387 
388  ATH_MSG_DEBUG ("station=" << station << ", plate= "<< plate << ", fiber=" << fiber << ", sign=" << sign << ", dep energy=" << ((*it)->GetEnergyDeposit()));
389 
390  m_E_fib[station-1][2*(plate-1)+(1-sign)/2][fiber-1] += ((*it)->GetEnergyDeposit());
391 
392  }
393 
394 }
395 
396 
397 
398 
399 void ALFA_PileUpTool::ALFA_MD_info(const ALFA_HitCollection* ALFA_HitCollection)
400 {
401 
402 // cleaning
403 
404  for ( int l = 0; l < 8; l++ )
405  {
406  for ( int i = 0; i < 20; i++ )
407  {
408  for ( int j = 0; j < 64; j++ )
409  {
410 
411  m_E_fib[l][i][j] = 0.;
412  m_fibres[l][i][j] = 0;
413 
414 
415  }
416  }
417  }
418 
419 // cleaning end
420 
421  int fiber, plate, sign, station;
422 
423  ALFA_HitConstIter it = ALFA_HitCollection->begin();
424  ALFA_HitConstIter itend = ALFA_HitCollection->end();
425 
426  for (; it != itend; ++it) {
427 
428  station = it->GetStationNumber();
429  plate = it->GetPlateNumber();
430  fiber = it->GetFiberNumber();
431  sign = it->GetSignFiber();
432 
433  ATH_MSG_DEBUG ("station=" << station << ", plate= "<< plate << ", fiber=" << fiber << ", sign=" << sign << ", dep energy=" << it->GetEnergyDeposit());
434 
435  m_E_fib[station-1][2*(plate-1)+(1-sign)/2][fiber-1] += it->GetEnergyDeposit();
436 
437 
438  }
439 }
440 
441 
442 
443 
444 void ALFA_PileUpTool::ALFA_OD_info(TimedHitCollection<ALFA_ODHit>& tODHitALFA)
445 {
446 
447 // cleaning
448 
449  for ( int l = 0; l < 8; l++ ){
450  for ( int i = 0; i < 2; i++ ){
451  for ( int j = 0; j < 3; j++ ){
452  for ( int k = 0; k < 30; k++ ){
453 
454  m_E_ODfib[l][i][j][k] = 0.;
455 
456  }
457  }
458  }
459  }
460 
461 // cleaning end
462 
463  int fiber, plate, sign, side, station;
464 
465 
466  TimedHitCollection<ALFA_ODHit> thpc = tODHitALFA;
467  TimedHitCollection<ALFA_ODHit>::const_iterator i, e, it;
468 
469  while (thpc.nextDetectorElement(i, e)) for (it = i; it != e; ++it) {
470 
471  station = (*it)->GetStationNumber();
472  side = (*it)->GetODSide();
473  plate = (*it)->GetPlateNumber();
474  fiber = (*it)->GetFiberNumber();
475  sign = (*it)->GetSignFiber();
476 
477 
478  ATH_MSG_DEBUG ("station=" << station << ", side=" << side << ", plate= "<< plate << ", fiber=" << fiber << ", sign=" << sign << ", dep energy=" << ((*it)->GetEnergyDeposit()));
479 
480  if (sign==0) m_E_ODfib[station-1][side-1][plate-1][fiber+15] += ((*it)->GetEnergyDeposit());
481  else m_E_ODfib[station-1][side-1][plate-1][fiber] += ((*it)->GetEnergyDeposit());
482 
483  }
484 
485 }
486 
487 
488 
489 
490 void ALFA_PileUpTool::ALFA_OD_info(const ALFA_ODHitCollection* ALFA_ODHitCollection)
491 {
492 
493  // cleaning
494 
495  for ( int l = 0; l < 8; l++ ){
496  for ( int i = 0; i < 2; i++ ){
497  for ( int j = 0; j < 3; j++ ){
498  for ( int k = 0; k < 30; k++ ){
499 
500  m_E_ODfib[l][i][j][k] = 0.;
501 
502  }
503  }
504  }
505  }
506 
507 // cleaning end
508 
509  int fiber, plate, sign, side, station;
510 
511 
512  ALFA_ODHitConstIter it = ALFA_ODHitCollection->begin();
513  ALFA_ODHitConstIter itend = ALFA_ODHitCollection->end();
514 
515  for (; it != itend; ++it) {
516 
517 
518  station = it->GetStationNumber();
519  side = it->GetODSide();
520  plate = it->GetPlateNumber();
521  fiber = it->GetFiberNumber();
522  sign = it->GetSignFiber();
523 
524 
525  ATH_MSG_DEBUG ("station=" << station << ", side=" << side << ", plate= "<< plate << ", fiber=" << fiber << ", sign=" << sign << ", dep energy=" << it->GetEnergyDeposit());
526 
527  if (sign==0) m_E_ODfib[station-1][side-1][plate-1][fiber+15] += it->GetEnergyDeposit();
528  else m_E_ODfib[station-1][side-1][plate-1][fiber] += it->GetEnergyDeposit();
529 
530  }
531 }
532 
533 
534 
535 
536 StatusCode ALFA_PileUpTool::fill_MD_DigitCollection(CLHEP::HepRandomEngine* rndEngine) {
537 
538  ATH_MSG_DEBUG(" ALFA_PileUpTool::fill_MD_DigitCollection()");
539 
540  double N_photo = 0.;
541  //double N_CTphoto = 0.;
542  //double sigma = 0.;
543  double amplitude = 0.;
544  //double amplitude_CT = 0.;
545 
546 
547  for ( int l = 0; l < 8; l++ )
548  {
549  for ( int i = 0; i < 20; i++ )
550  {
551  for ( int j = 0; j < 64; j++ )
552  {
553 
554  N_photo = CLHEP::RandPoisson::shoot(rndEngine,m_E_fib[l][i][j]*m_meanN_photo/m_meanE_dep);
555 
556  //sigma = m_sigma1 * sqrt(N_photo);
557  //double noise_1 = sigma * CLHEP::RandGaussZiggurat::shoot (rndEngine, m_mean, m_stdDev);
558  //double noise_2 = m_sigma0 * CLHEP::RandGaussZiggurat::shoot (rndEngine, m_mean, m_stdDev);
559 
560  //amplitude = N_photo
561  // + noise_1
562  // + noise_2;
563 
564  amplitude = CLHEP::RandGaussZiggurat::shoot(rndEngine, N_photo, sqrt(pow(m_sigma0,2)+N_photo*pow(m_sigma1,2)));
565 
566  //N_photo_CT = CLHEP::RandPoisson::shoot(rndEngine,0.08*m_E_fib[l][i][j]*m_meanN_photo/m_meanE_dep);
567  //amplitude_CT = CLHEP::RandGaussQ::shoot (rndEngine, N_photo_CT, sqrt(pow(m_sigma0,2)+N_photo_CT*pow(m_sigma1,2)));
568 
569 
570  if (amplitude >= m_AmplitudeCut )
571  {
572  ATH_MSG_DEBUG(" ALFA_Digitization::fillDigitCollection, amplitude " << amplitude);
573  ATH_MSG_DEBUG(" station = " << l << ", plate= " << i << ", fiber=" << j );
574  m_digitCollection->push_back(new ALFA_Digit(l,i,j));
575  m_fibres[l][i][j] = m_fibres[l][i][j] + 1;
576  }
577 
578  }
579  }
580  }
581 
582  double rand_fib;
583 
584  for ( int l = 0; l < 8; l++ )
585  {
586  for ( int i = 0; i < 20; i++ )
587  {
588  for ( int j = 0; j < 64; j++ )
589  {
590  if ( m_fibres[l][i][j] > 0 )
591  {
592  for (int f = j+1; f < 64; f++)
593  {
594  rand_fib = CLHEP::RandFlat::shoot(rndEngine,0.,1.);
595 
596  if ( m_fibres[l][i][f] == 0)
597  {
598  if (rand_fib <= m_cross_talk[l][63+f-j])
599  {
600  m_digitCollection->push_back(new ALFA_Digit(l,i,f));
601  }
602  }
603  }
604 
605  for (int f = j-1; f > -1; f--)
606  {
607  rand_fib = CLHEP::RandFlat::shoot(rndEngine,0.,1.);
608 
609  if ( m_fibres[l][i][f] == 0)
610  {
611  if (rand_fib <= m_cross_talk[l][63-(j-f)])
612  {
613  m_digitCollection->push_back(new ALFA_Digit(l,i,f));
614  }
615  }
616  }
617  }
618  }
619  }
620  }
621  return StatusCode::SUCCESS;
622 
623 }
624 
625 
626 
627 
628 StatusCode ALFA_PileUpTool::fill_OD_DigitCollection(CLHEP::HepRandomEngine* rndEngine) {
629 
630  ATH_MSG_DEBUG(" ALFA_PileUpTool::fill_OD_DigitCollection ");
631 
632  double N_photo = 0.;
633  double sigma = 0.;
634  double amplitude = 0.;
635 
636  for ( int l = 0; l < 8; l++ )
637  {
638  for ( int i = 0; i < 2; i++ )
639  {
640  for ( int j = 0; j < 3; j++ )
641  {
642  for (int k = 0; k < 30; k++)
643  {
644 
645  N_photo = CLHEP::RandPoisson::shoot(rndEngine,m_E_ODfib[l][i][j][k]*m_meanN_photo/m_meanE_dep);
646  sigma = m_sigma1 * sqrt(N_photo);
647 
648  double noise_1 = sigma * CLHEP::RandGaussZiggurat::shoot (rndEngine, m_mean, m_stdDev);
649  double noise_2 = m_sigma0 * CLHEP::RandGaussZiggurat::shoot (rndEngine, m_mean, m_stdDev);
650 
651 
652  amplitude = N_photo
653  + noise_1
654  + noise_2;
655 
656 
657  if (amplitude >= m_AmplitudeCut )
658  {
659  ATH_MSG_DEBUG(" ALFA_Digitization::fill_OD_DigitCollection, amplitude " << amplitude);
660  ATH_MSG_DEBUG(" station = " << l << ", side = " << i << ", plate = " << j << ", fiber = " << k);
661 
662  m_ODdigitCollection->push_back(new ALFA_ODDigit(l,i,j,k));
663  }
664 
665 
666  }
667  }
668  }
669  }
670 
671  return StatusCode::SUCCESS;
672 }
673 
674 
675  StatusCode ALFA_PileUpTool::XTalk()
676  {
677  for (unsigned int j=0; j<8; j++){
678 
679  m_s.str("");
680 
681  m_filename = "Xtalk_station";
682  m_s << j+1;
683  m_filename += m_s.str();
684  m_filename += ".txt";
685 
686  ATH_MSG_DEBUG("file name " << m_filename.c_str() );
687 
688  std::string filePath = PathResolver::find_file(m_filename,"DATAPATH", PathResolver::RecursiveSearch);
689 
690  if(filePath.length() == 0)
691  {
692  ATH_MSG_FATAL(" XTalk file " << m_filename.c_str() << " not found in Datapath");
693  throw std::runtime_error("FATAL: mapping MD maroc-mapmt not found in Datapath.");
694  }
695 
696  else
697  {
698  ATH_MSG_DEBUG("the XTALK file \"" << m_filename.c_str() << "\" found in Datapath");
699  ATH_MSG_DEBUG("filePath = " << filePath.c_str() );
700  }
701 
702  m_fXTalk.open(filePath.c_str());
703 
704  if (m_fXTalk.is_open())
705  {
706  for (unsigned int i=0;i<127;i++)
707  {
708  m_fXTalk >> m_cross_talk[j][i];
709 
710  ATH_MSG_DEBUG("cross_talk[" << j << "][" << i << "] = " << m_cross_talk[j][i]);
711 
712  }
713  m_fXTalk.close();
714  }
715  else
716  {
717  ATH_MSG_ERROR("the file " << m_filename.c_str() << " was not open");
718  return StatusCode::FAILURE;
719  }
720  }
721  return StatusCode::SUCCESS;
722  }
xAOD::iterator
JetConstituentVector::iterator iterator
Definition: JetConstituentVector.cxx:68
ALFA_PileUpTool::m_fibres
int m_fibres[8][20][64]
Definition: ALFA_PileUpTool.h:163
ATHRNG::RNGWrapper::setSeed
void setSeed(const std::string &algName, const EventContext &ctx)
Set the random seed using a string (e.g.
Definition: RNGWrapper.h:169
ALFA_HitCollection
AtlasHitsVector< ALFA_Hit > ALFA_HitCollection
Definition: ALFA_HitCollection.h:13
ALFA_PileUpTool::ALFA_MD_info
void ALFA_MD_info(const ALFA_HitCollection *)
Definition: ALFA_PileUpTool.cxx:399
python.CaloRecoConfig.f
f
Definition: CaloRecoConfig.py:127
ALFA_ODDigitCollection.h
ATH_MSG_FATAL
#define ATH_MSG_FATAL(x)
Definition: AthMsgStreamMacros.h:34
pdg_comparison.sigma
sigma
Definition: pdg_comparison.py:324
GenEvent.h
PathResolver::RecursiveSearch
@ RecursiveSearch
Definition: PathResolver.h:28
ALFA_ODHit.h
ALFA_PileUpTool::processBunchXing
virtual StatusCode processBunchXing(int bunchXing, SubEventIterator bSubEvents, SubEventIterator eSubEvents) override final
called for each active bunch-crossing to process current SubEvents bunchXing is in ns
Definition: ALFA_PileUpTool.cxx:243
ALFA_PileUpTool::m_meanN_photo
double m_meanN_photo
Definition: ALFA_PileUpTool.h:146
PathResolver::find_file
static std::string find_file(const std::string &logical_file_name, const std::string &search_path, SearchType search_type=LocalSearch)
Definition: PathResolver.cxx:251
AthMsgStreamMacros.h
ALFA_PileUpTool::m_key_ODDigitCollection
std::string m_key_ODDigitCollection
Definition: ALFA_PileUpTool.h:127
GenVertex.h
conifer::pow
constexpr int pow(int x)
Definition: conifer.h:20
AtlasHitsVector
Definition: AtlasHitsVector.h:33
ALFA_ODHitCollection.h
skel.it
it
Definition: skel.GENtoEVGEN.py:423
ALFA_PileUpTool::ALFA_OD_info
void ALFA_OD_info(const ALFA_ODHitCollection *)
Definition: ALFA_PileUpTool.cxx:490
UploadAMITag.l
list l
Definition: UploadAMITag.larcaf.py:158
ALFA_PileUpTool::m_filename
std::string m_filename
Definition: ALFA_PileUpTool.h:167
ALFA_HitConstIter
AtlasHitsVector< ALFA_Hit >::const_iterator ALFA_HitConstIter
Definition: ALFA_HitCollection.h:15
ATH_MSG_VERBOSE
#define ATH_MSG_VERBOSE(x)
Definition: AthMsgStreamMacros.h:28
AtlasHitsVector::begin
const_iterator begin() const
Definition: AtlasHitsVector.h:131
ALFA_PileUpTool::m_AmplitudeCut
double m_AmplitudeCut
Definition: ALFA_PileUpTool.h:149
TimedHitCollection::nextDetectorElement
bool nextDetectorElement(const_iterator &b, const_iterator &e)
sets an iterator range with the hits of current detector element returns a bool when done
GenParticle.h
AtlasHitsVector::const_iterator
CONT::const_iterator const_iterator
Definition: AtlasHitsVector.h:43
ALFA_PileUpTool.h
AthenaPoolTestRead.sc
sc
Definition: AthenaPoolTestRead.py:27
AtlasHitsVector::clear
void clear()
Definition: AtlasHitsVector.h:153
ALFA_PileUpTool::initialize
virtual StatusCode initialize() override final
Definition: ALFA_PileUpTool.cxx:93
PileUpMergeSvc::TimedList::type
std::list< value_t > type
type of the collection of timed data object
Definition: PileUpMergeSvc.h:75
StoreGateSvc::retrieve
StatusCode retrieve(const T *&ptr) const
Retrieve the default object into a const T*.
TRT::Hit::side
@ side
Definition: HitInfo.h:83
ALFA_PileUpTool::m_mergedALFA_HitList
ALFA_HitCollection * m_mergedALFA_HitList
Definition: ALFA_PileUpTool.h:159
ALFA_PileUpTool::m_randomStreamName
Gaudi::Property< std::string > m_randomStreamName
Definition: ALFA_PileUpTool.h:113
StoreGateSvc
The Athena Transient Store API.
Definition: StoreGateSvc.h:128
ALFA_PileUpTool::m_ODdigitCollection
ALFA_ODDigitCollection * m_ODdigitCollection
Definition: ALFA_PileUpTool.h:132
ALFA_Digit
Definition: ALFA_Digit.h:8
SimpleVector.h
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition: AthMsgStreamMacros.h:33
ALFA_PileUpTool::recordCollection
StatusCode recordCollection(ServiceHandle< StoreGateSvc > &evtStore, const std::string &key_digitCnt)
Definition: ALFA_PileUpTool.cxx:329
TimedHitCollection::insert
void insert(const PileUpTimeEventIndex &timeEventIndex, const AtlasHitsVector< HIT > *inputCollection)
McEventCollection.h
lumiFormat.i
int i
Definition: lumiFormat.py:92
ALFA_PileUpTool::m_meanE_dep
double m_meanE_dep
Definition: ALFA_PileUpTool.h:145
EL::StatusCode
::StatusCode StatusCode
StatusCode definition for legacy code.
Definition: PhysicsAnalysis/D3PDTools/EventLoop/EventLoop/StatusCode.h:22
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition: AthMsgStreamMacros.h:29
ALFA_PileUpTool::fill_OD_DigitCollection
StatusCode fill_OD_DigitCollection(CLHEP::HepRandomEngine *)
Definition: ALFA_PileUpTool.cxx:628
ALFA_PileUpTool::m_SimODHitCollectionName
std::string m_SimODHitCollectionName
Definition: ALFA_PileUpTool.h:125
ALFA_PileUpTool::m_E_ODfib
double m_E_ODfib[8][2][3][30]
Definition: ALFA_PileUpTool.h:116
test_pyathena.parent
parent
Definition: test_pyathena.py:15
sign
int sign(int a)
Definition: TRT_StrawNeighbourSvc.h:127
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
ALFA_ODHitCollection
AtlasHitsVector< ALFA_ODHit > ALFA_ODHitCollection
Definition: ALFA_ODHitCollection.h:13
ALFA_DigitCollection.h
ALFA_PileUpTool::m_E_fib
double m_E_fib[8][20][64]
Definition: ALFA_PileUpTool.h:115
PileUpToolBase
Definition: PileUpToolBase.h:18
ALFA_PileUpTool::processAllSubEvents
virtual StatusCode processAllSubEvents(const EventContext &ctx) override final
Definition: ALFA_PileUpTool.cxx:116
ALFA_PileUpTool::m_s
std::stringstream m_s
Definition: ALFA_PileUpTool.h:166
ALFA_PileUpTool::recordODCollection
StatusCode recordODCollection(ServiceHandle< StoreGateSvc > &evtStore, const std::string &key_ODdigitCnt)
Definition: ALFA_PileUpTool.cxx:341
ALFA_PileUpTool::m_mergedALFA_ODHitList
ALFA_ODHitCollection * m_mergedALFA_ODHitList
Definition: ALFA_PileUpTool.h:160
ALFA_PileUpTool::m_sigma0
double m_sigma0
Definition: ALFA_PileUpTool.h:143
PathResolver.h
ALFA_DigitCollection
Definition: ALFA_DigitCollection.h:15
ALFA_PileUpTool::XTalk
StatusCode XTalk()
Definition: ALFA_PileUpTool.cxx:675
ALFA_PileUpTool::m_sigma1
double m_sigma1
Definition: ALFA_PileUpTool.h:144
ALFA_Hit.h
hancool.filePath
string filePath
Definition: hancool.py:28
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:195
ATHRNG::RNGWrapper
A wrapper class for event-slot-local random engines.
Definition: RNGWrapper.h:56
ALFA_PileUpTool::m_SimHitCollectionName
std::string m_SimHitCollectionName
Definition: ALFA_PileUpTool.h:124
ALFA_ODHitConstIter
AtlasHitsVector< ALFA_ODHit >::const_iterator ALFA_ODHitConstIter
Definition: ALFA_ODHitCollection.h:15
DataVector::push_back
value_type push_back(value_type pElem)
Add an element to the end of the collection.
ALFA_PileUpTool::finalize
virtual StatusCode finalize() override final
Definition: ALFA_PileUpTool.cxx:325
ALFA_PileUpTool::m_cross_talk
double m_cross_talk[8][127]
Definition: ALFA_PileUpTool.h:162
ATHRNG::RNGWrapper::getEngine
CLHEP::HepRandomEngine * getEngine(const EventContext &ctx) const
Retrieve the random engine corresponding to the provided EventContext.
Definition: RNGWrapper.h:134
RNGWrapper.h
ALFA_PileUpTool::m_mergeSvc
ServiceHandle< PileUpMergeSvc > m_mergeSvc
Definition: ALFA_PileUpTool.h:111
ALFA_PileUpTool::m_key_DigitCollection
std::string m_key_DigitCollection
Definition: ALFA_PileUpTool.h:126
DiTauMassTools::MaxHistStrategyV2::e
e
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:26
ALFA_PileUpTool::m_fXTalk
std::ifstream m_fXTalk
Definition: ALFA_PileUpTool.h:165
ALFA_PileUpTool::prepareEvent
virtual StatusCode prepareEvent(const EventContext &ctx, const unsigned int nInputEvents) override final
code taken from ZDC; author (highly) probably John Chapman called before the subevts loop.
Definition: ALFA_PileUpTool.cxx:222
ALFA_HitCollection.h
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition: AthMsgStreamMacros.h:32
AtlasHitsVector::end
const_iterator end() const
Definition: AtlasHitsVector.h:134
ALFA_PileUpTool::ALFA_PileUpTool
ALFA_PileUpTool(const std::string &type, const std::string &name, const IInterface *parent)
Definition: ALFA_PileUpTool.cxx:49
python.CaloScaleNoiseConfig.type
type
Definition: CaloScaleNoiseConfig.py:78
TimedHitCollection::const_iterator
TimedVector::const_iterator const_iterator
Definition: TimedHitCollection.h:20
ALFA_PileUpTool::m_digitCollection
ALFA_DigitCollection * m_digitCollection
Definition: ALFA_PileUpTool.h:130
declareProperty
#define declareProperty(n, p, h)
Definition: BaseFakeBkgTool.cxx:15
AtlasHitsVector::size
size_type size() const
Definition: AtlasHitsVector.h:143
SubEventIterator
std::vector< xAOD::EventInfo::SubEvent >::const_iterator SubEventIterator
Definition: IPileUpTool.h:22
ALFA_PileUpTool::m_mean
double m_mean
Definition: ALFA_PileUpTool.h:147
ALFA_PileUpTool::mergeEvent
virtual StatusCode mergeEvent(const EventContext &ctx) override final
return false if not interested in certain xing times (in ns) implemented by default in PileUpToolBase...
Definition: ALFA_PileUpTool.cxx:298
ALFA_ODDigit
Definition: ALFA_ODDigit.h:8
ALFA_PileUpTool::fill_MD_DigitCollection
StatusCode fill_MD_DigitCollection(CLHEP::HepRandomEngine *)
Definition: ALFA_PileUpTool.cxx:536
ALFA_PileUpTool::m_randomSvc
ServiceHandle< IAthRNGSvc > m_randomSvc
Definition: ALFA_PileUpTool.h:112
TimedHitCollection
Definition: TimedHitCollection.h:15
AtlasHitsVector::push_back
void push_back(const T &t)
Definition: AtlasHitsVector.h:145
fitman.k
k
Definition: fitman.py:528
ALFA_PileUpTool::m_stdDev
double m_stdDev
Definition: ALFA_PileUpTool.h:148
ServiceHandle< StoreGateSvc >
ALFA_ODDigitCollection
Definition: ALFA_ODDigitCollection.h:15
Generated on Fri Jul 5 2024 21:10:10 for ATLAS Offline Software by   doxygen 1.8.18