gambit is hosted by Hepforge, IPPP Durham
GAMBIT  v1.5.0-2191-ga4742ac
a Global And Modular Bsm Inference Tool
generateEventPy8Collider.hpp
Go to the documentation of this file.
1 // GAMBIT: Global and Modular BSM Inference Tool
2 // *********************************************
43 
46 
47 // #define COLLIDERBIT_DEBUG
48 #define DEBUG_PREFIX "DEBUG: OMP thread " << omp_get_thread_num() << ": "
49 
50 namespace Gambit
51 {
52 
53  namespace ColliderBit
54  {
55 
57  #ifndef EXCLUDE_HEPMC
58  template<typename PythiaT, typename hepmc_writerT>
59  void dropHepMCEventPy8Collider(const PythiaT* Pythia, const safe_ptr<Options>& runOptions)
60  {
61  // Write event to HepMC file
62  static const bool drop_HepMC2_file = runOptions->getValueOrDef<bool>(false, "drop_HepMC2_file");
63  static const bool drop_HepMC3_file = runOptions->getValueOrDef<bool>(false, "drop_HepMC3_file");
64  if (drop_HepMC2_file or drop_HepMC3_file)
65  {
66  thread_local hepmc_writerT hepmc_writer;
67  thread_local bool first = true;
68 
69  if (first)
70  {
71  str filename = "GAMBIT_collider_events.omp_thread_";
72  filename += std::to_string(omp_get_thread_num());
73  filename += ".hepmc";
74  hepmc_writer.init(filename, drop_HepMC2_file, drop_HepMC3_file);
75  first = false;
76  }
77 
78  if(drop_HepMC2_file)
79  hepmc_writer.write_event_HepMC2(const_cast<PythiaT*>(Pythia));
80  if(drop_HepMC3_file)
81  hepmc_writer.write_event_HepMC3(const_cast<PythiaT*>(Pythia));
82 
83  }
84  }
85  #endif
86 
88  template<typename PythiaT, typename EventT, typename hepmc_writerT>
89  void generateEventPy8Collider(HEPUtils::Event& event,
90  const MCLoopInfo& RunMC,
91  const Py8Collider<PythiaT,EventT,hepmc_writerT>& HardScatteringSim,
92  const EventWeighterFunctionType& EventWeighterFunction,
93  const int iteration,
94  void(*wrapup)(),
95  const safe_ptr<Options>& runOptions)
96  {
97  static int nFailedEvents;
98  thread_local EventT pythia_event;
99  static const double jet_pt_min = runOptions->getValueOrDef<double>(10.0, "jet_pt_min");
100 
101  // If the event loop has not been entered yet, reset the counter for the number of failed events
102  if (iteration == BASE_INIT)
103  {
104  // Although BASE_INIT should never be called in parallel, we do this in omp_critical just in case.
105  #pragma omp critical (pythia_event_failure)
106  {
107  nFailedEvents = 0;
108  }
109  return;
110  }
111 
112  // If in any other special iteration, do nothing
113  if (iteration < BASE_INIT) return;
114 
115  // Reset the Pythia and HEPUtils events
116  pythia_event.clear();
117  event.clear();
118 
119  // Attempt (possibly repeatedly) to generate an event
120  while(nFailedEvents <= RunMC.current_maxFailedEvents())
121  {
122  try
123  {
124  #ifdef COLLIDERBIT_DEBUG
125  cerr << DEBUG_PREFIX << "Will now call HardScatteringSim.nextEvent(pythia_event)..." << endl;
126  #endif
127 
128  HardScatteringSim.nextEvent(pythia_event);
129  break;
130  }
132  {
133  #ifdef COLLIDERBIT_DEBUG
134  cerr << DEBUG_PREFIX << "Py8Collider::EventGenerationError caught in generateEventPy8Collider. Check the ColliderBit log for event details." << endl;
135  #endif
136  #pragma omp critical (pythia_event_failure)
137  {
138  // Update global counter
139  nFailedEvents += 1;
140  // Store Pythia event record in the logs
141  std::stringstream ss;
142  pythia_event.list(ss, 1);
143  logger() << LogTags::debug << "Py8Collider::EventGenerationError error caught in generateEventPy8Collider. Pythia record for event that failed:\n" << ss.str() << EOM;
144  }
145  }
146  }
147 
148  // Wrap up event loop if too many events fail.
149  if(nFailedEvents > RunMC.current_maxFailedEvents())
150  {
151  // Tell the MCLoopInfo instance that we have exceeded maxFailedEvents
154  {
155  piped_invalid_point.request("exceeded maxFailedEvents");
156  }
157  wrapup();
158  return;
159  }
160 
161  #ifndef EXCLUDE_HEPMC
162  dropHepMCEventPy8Collider<PythiaT,hepmc_writerT>(HardScatteringSim.pythia(), runOptions);
163  #endif
164 
165 
166  // Attempt to convert the Pythia event to a HEPUtils event
167  try
168  {
169  if (HardScatteringSim.partonOnly)
170  convertPartonEvent(pythia_event, event, HardScatteringSim.antiktR, jet_pt_min);
171  else
172  convertParticleEvent(pythia_event, event, HardScatteringSim.antiktR, jet_pt_min);
173  }
174  // No good.
175  catch (Gambit::exception& e)
176  {
177  #ifdef COLLIDERBIT_DEBUG
178  cerr << DEBUG_PREFIX << "Gambit::exception caught during event conversion in generateEventPy8Collider. Check the ColliderBit log for details." << endl;
179  #endif
180 
181  #pragma omp critical (event_conversion_error)
182  {
183  // Store Pythia event record in the logs
184  std::stringstream ss;
185  pythia_event.list(ss, 1);
186  logger() << LogTags::debug << "Gambit::exception caught in generateEventPy8Collider. Pythia record for event that failed:\n" << ss.str() << EOM;
187  }
188 
189  str errmsg = "Bad point: generateEventPy8Collider caught the following runtime error: ";
190  errmsg += e.what();
192  wrapup();
193  return;
194  }
195 
196  // Assign weight to event
197  EventWeighterFunction(event, &HardScatteringSim);
198  }
199 
201  #define GET_PYTHIA_EVENT(NAME) \
202  void NAME(HEPUtils::Event& result) \
203  { \
204  using namespace Pipes::NAME; \
205  generateEventPy8Collider(result, *Dep::RunMC, \
206  *Dep::HardScatteringSim, *Dep::EventWeighterFunction, \
207  *Loop::iteration, Loop::wrapup,runOptions); \
208  }
209 
210  }
211 
212 }
void generateEventPy8Collider(HEPUtils::Event &event, const MCLoopInfo &RunMC, const Py8Collider< PythiaT, EventT, hepmc_writerT > &HardScatteringSim, const EventWeighterFunctionType &EventWeighterFunction, const int iteration, void(*wrapup)(), const safe_ptr< Options > &runOptions)
Generate a hard scattering event with Pythia.
void report_exceeded_maxFailedEvents() const
Set exceeded_maxFailedEvents = true and decrement event counter by 1.
Definition: MCLoopInfo.cpp:29
A safe pointer that throws an informative error if you try to dereference it when nullified...
Definition: util_types.hpp:175
const int & current_maxFailedEvents() const
Get maximum allowable number of failed events before MC loop is terminated for the current collider...
Definition: MCLoopInfo.cpp:129
bool partonOnly
Flag indicating if events from this collider should be processed as parton-only or full events...
void nextEvent(EventT &event) const
Event generation for any Pythia interface to Gambit.
double antiktR
The jet radius used for the anti-kt jet clustering.
Piped_invalid_point piped_invalid_point
Global instance of piped invalid point class.
Definition: exceptions.cpp:544
Declarations common to all ColliderBit event loop functions.
void convertPartonEvent(const EventT &pevt, HEPUtils::Event &result, double antiktR, double jet_pt_min)
Convert a partonic (no hadrons) EventT into an unsmeared HEPUtils::Event.
const Logging::endofmessage EOM
Explicit const instance of the end of message struct in Gambit namespace.
Definition: logger.hpp:100
EXPORT_SYMBOLS Logging::LogMaster & logger()
Function to retrieve a reference to the Gambit global log object.
Definition: logger.cpp:95
const bool & current_invalidate_failed_points() const
Get invalidate_failed_points bool for the current collider.
Definition: MCLoopInfo.cpp:132
void dropHepMCEventPy8Collider(const PythiaT *Pythia, const safe_ptr< Options > &runOptions)
Drop a HepMC file for the event.
std::string str
Shorthand for a standard string.
Definition: Analysis.hpp:35
#define DEBUG_PREFIX
An exception for when Pythia fails to generate events.
Definition: Py8Collider.hpp:71
const PythiaT * pythia() const
Get the Pythia instance.
Definition: Py8Collider.hpp:49
void convertParticleEvent(const EventT &pevt, HEPUtils::Event &result, double antiktR, double jet_pt_min)
Convert a hadron-level EventT into an unsmeared HEPUtils::Event.
A specializable, recyclable class interfacing ColliderBit and Pythia.
Definition: Py8Collider.hpp:36
void request(std::string message)
Request an exception.
Definition: exceptions.cpp:471
std::function< void(HEPUtils_Event &, const BaseCollider *)> EventWeighterFunctionType
TODO: see if we can use this one:
Definition: Analysis.hpp:33
Helper functions for converting between Pythia8 events and other event types.
Container for event loop status data and settings.
Definition: MCLoopInfo.hpp:31