gambit is hosted by Hepforge, IPPP Durham
GAMBIT  v1.5.0-2191-ga4742ac
a Global And Modular Bsm Inference Tool
Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT > Class Template Reference

A specializable, recyclable class interfacing ColliderBit and Pythia. More...

#include <Py8Collider.hpp>

Inheritance diagram for Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >:
Collaboration diagram for Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >:

Classes

class  EventGenerationError
 An exception for when Pythia fails to generate events. More...
 
class  InitializationError
 An exception for when Pythia fails to initialize. More...
 

Public Member Functions

const PythiaT * pythia () const
 Get the Pythia instance. More...
 
bool SetupMatchingUserHook ()
 
Construction, Destruction, and Recycling:
 Py8Collider ()
 
 ~Py8Collider ()
 
void clear ()
 Reset this instance for reuse, avoiding the need for "new" or "delete". More...
 
(Re-)Initialization functions
void addToSettings (const std::string &command)
 Add a command to the list of settings used by "init". More...
 
void banner (const std::string pythiaDocPath)
 Create a useless Pythia instance just to print the banner. More...
 
void init ()
 Initialize with no settings (error): override version. More...
 
void init (const std::vector< std::string > &externalSettings)
 Initialize from some external settings: override version. More...
 
void init (const std::string pythiaDocPath, const std::vector< std::string > &externalSettings, const SLHAea::Coll *slhaea=nullptr, std::ostream &os=std::cout)
 Initialize from some external settings. More...
 
void init_user_model (const std::string pythiaDocPath, const std::vector< std::string > &externalSettings, const SLHAea::Coll *slhaea=nullptr, std::ostream &os=std::cout)
 Initialize from some external settings. More...
 
void init (const std::string pythiaDocPath, const std::vector< std::string > &externalSettings, std::ostream &os)
 Initialize from some external settings, assuming no given SLHAea instance. More...
 
void init_user_model (const std::string pythiaDocPath, const std::vector< std::string > &externalSettings, std::ostream &os)
 Initialize from some external settings, assuming no given SLHAea instance. More...
 
Event generation and cross section functions
void nextEvent (EventT &event) const
 Event generation for any Pythia interface to Gambit. More...
 
double xsec_fb () const
 Report the total or process-specific cross section (in fb or pb). More...
 
double xsec_fb (int process_code) const
 
double xsec_pb () const
 
double xsec_pb (int process_code) const
 
double xsecErr_fb () const
 Report the uncertainty in the total or process-specific cross section (in fb or pb). More...
 
double xsecErr_fb (int process_code) const
 
double xsecErr_pb () const
 
double xsecErr_pb (int process_code) const
 
int process_code () const
 Report an integer process code for the last generated event. More...
 
std::vector< intall_active_process_codes () const
 Report the list of all active process codes. More...
 
- Public Member Functions inherited from Gambit::ColliderBit::BaseCollider
 BaseCollider ()
 Constructor. More...
 
virtual ~BaseCollider ()
 Destructor. More...
 

Protected Attributes

PythiaT * _pythiaInstance
 
PythiaT * _pythiaBase
 
std::vector< std::string > _pythiaSettings
 

Additional Inherited Members

- Public Attributes inherited from Gambit::ColliderBit::BaseCollider
bool partonOnly
 Flag indicating if events from this collider should be processed as parton-only or full events. More...
 
double antiktR
 The jet radius used for the anti-kt jet clustering. More...
 

Detailed Description

template<typename PythiaT, typename EventT, typename hepmc_writerT>
class Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >

A specializable, recyclable class interfacing ColliderBit and Pythia.

Definition at line 36 of file Py8Collider.hpp.

Constructor & Destructor Documentation

◆ Py8Collider()

template<typename PythiaT, typename EventT, typename hepmc_writerT>
Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::Py8Collider ( )
inline

Definition at line 85 of file Py8Collider.hpp.

◆ ~Py8Collider()

template<typename PythiaT, typename EventT, typename hepmc_writerT>
Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::~Py8Collider ( )
inline

Member Function Documentation

◆ addToSettings()

template<typename PythiaT, typename EventT, typename hepmc_writerT>
void Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::addToSettings ( const std::string &  command)
inline

Add a command to the list of settings used by "init".

Definition at line 111 of file Py8Collider.hpp.

111 { _pythiaSettings.push_back(command); }
std::vector< std::string > _pythiaSettings
Definition: Py8Collider.hpp:43

◆ all_active_process_codes()

template<typename PythiaT, typename EventT, typename hepmc_writerT>
std::vector<int> Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::all_active_process_codes ( ) const
inlinevirtual

Report the list of all active process codes.

Implements Gambit::ColliderBit::BaseCollider.

Definition at line 237 of file Py8Collider.hpp.

237 { return _pythiaInstance->info.codesHard(); }

◆ banner()

template<typename PythiaT, typename EventT, typename hepmc_writerT>
void Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::banner ( const std::string  pythiaDocPath)
inline

Create a useless Pythia instance just to print the banner.

Definition at line 113 of file Py8Collider.hpp.

Referenced by Gambit::ColliderBit::getPy8Collider().

113 { PythiaT myPythia(pythiaDocPath); }
Here is the caller graph for this function:

◆ clear()

template<typename PythiaT, typename EventT, typename hepmc_writerT>
void Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::clear ( )
inlinevirtual

Reset this instance for reuse, avoiding the need for "new" or "delete".

Reimplemented from Gambit::ColliderBit::BaseCollider.

Definition at line 94 of file Py8Collider.hpp.

References Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::_pythiaInstance.

Referenced by Gambit::ColliderBit::getPy8Collider().

95  {
96  _pythiaSettings.clear();
97  if (_pythiaInstance)
98  {
99  delete _pythiaInstance;
100  _pythiaInstance=nullptr;
101  }
102  }
std::vector< std::string > _pythiaSettings
Definition: Py8Collider.hpp:43
Here is the caller graph for this function:

◆ init() [1/4]

template<typename PythiaT, typename EventT, typename hepmc_writerT>
void Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::init ( )
inlinevirtual

Initialize with no settings (error): override version.

Reimplemented from Gambit::ColliderBit::BaseCollider.

Definition at line 115 of file Py8Collider.hpp.

Referenced by Gambit::ColliderBit::getPy8Collider(), and Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::init().

115 { std::cout<<"No settings given to Pythia!\n\n"; throw InitializationError(); }
Here is the caller graph for this function:

◆ init() [2/4]

template<typename PythiaT, typename EventT, typename hepmc_writerT>
void Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::init ( const std::vector< std::string > &  externalSettings)
inlinevirtual

Initialize from some external settings: override version.

Note
A string denoting the path to Pythia's xmldoc directory is
assumed to be at the end of the settings vector:

Reimplemented from Gambit::ColliderBit::BaseCollider.

Definition at line 120 of file Py8Collider.hpp.

References Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::init().

121  {
122  std::string docPath = externalSettings.back();
123  std::vector<std::string> settings(externalSettings);
124  settings.pop_back();
125  init(docPath, settings);
126  }
void init()
Initialize with no settings (error): override version.
Here is the call graph for this function:

◆ init() [3/4]

template<typename PythiaT, typename EventT, typename hepmc_writerT>
void Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::init ( const std::string  pythiaDocPath,
const std::vector< std::string > &  externalSettings,
const SLHAea::Coll *  slhaea = nullptr,
std::ostream &  os = std::cout 
)
inline

Initialize from some external settings.

Note
This override is most commonly used in ColliderBit.

Definition at line 130 of file Py8Collider.hpp.

References Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::_pythiaInstance.

133  {
134  // Settings acquired externally (ex from a gambit yaml file)
135  for(const auto command : externalSettings) _pythiaSettings.push_back(command);
136 
137  if (!_pythiaBase)
138  {
139  _pythiaBase = new PythiaT(pythiaDocPath, false);
140  }
141 
142  // Pass all settings to _pythiaBase
143  for(const auto command : _pythiaSettings) _pythiaBase->readString(command);
144 
145  // Create new _pythiaInstance from _pythiaBase
146  if (_pythiaInstance) delete _pythiaInstance;
147  _pythiaInstance = new PythiaT(_pythiaBase->particleData, _pythiaBase->settings);
148 
149  // Send along the SLHAea::Coll pointer, if it exists
150  if (slhaea) _pythiaInstance->slhaInterface.slha.setSLHAea(slhaea);
151 
152  // Read command again to get SM decay table change from yaml file
153  for(const auto command : _pythiaSettings)
154  {
155  _pythiaInstance->readString(command);
156  }
157 
158  if (!_pythiaInstance->init(os)) throw InitializationError();
159  }
std::vector< std::string > _pythiaSettings
Definition: Py8Collider.hpp:43

◆ init() [4/4]

template<typename PythiaT, typename EventT, typename hepmc_writerT>
void Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::init ( const std::string  pythiaDocPath,
const std::vector< std::string > &  externalSettings,
std::ostream &  os 
)
inline

Initialize from some external settings, assuming no given SLHAea instance.

Definition at line 190 of file Py8Collider.hpp.

References Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::init().

192  {
193  init(pythiaDocPath, externalSettings, nullptr, os);
194  }
void init()
Initialize with no settings (error): override version.
Here is the call graph for this function:

◆ init_user_model() [1/2]

template<typename PythiaT, typename EventT, typename hepmc_writerT>
void Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::init_user_model ( const std::string  pythiaDocPath,
const std::vector< std::string > &  externalSettings,
const SLHAea::Coll *  slhaea = nullptr,
std::ostream &  os = std::cout 
)
inline

Initialize from some external settings.

Special version of the init function for user defined models Needs to directly construct the new matrix elements (rather than use flags)

Definition at line 164 of file Py8Collider.hpp.

References Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::_pythiaInstance.

Referenced by Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::init_user_model().

167  {
168  // Settings acquired externally (for example, from a gambit yaml file)
169  for(const auto command : externalSettings) _pythiaSettings.push_back(command);
170 
171  if (!_pythiaBase)
172  {
173  _pythiaBase = new PythiaT(pythiaDocPath, false);
174  }
175 
176  // Pass all settings to _pythiaBase
177  for(const auto command : _pythiaSettings) _pythiaBase->readString(command);
178 
179  // Create new _pythiaInstance from _pythiaBase
180  if (_pythiaInstance) delete _pythiaInstance;
181  _pythiaInstance = new PythiaT(_pythiaBase->particleData, _pythiaBase->settings);
182 
183  // Send along the SLHAea::Coll pointer, if it exists
184  if (slhaea) _pythiaInstance->slhaInterface.slha.setSLHAea(slhaea);
185 
186  if (!_pythiaInstance->init(os)) throw InitializationError();
187  }
std::vector< std::string > _pythiaSettings
Definition: Py8Collider.hpp:43
Here is the caller graph for this function:

◆ init_user_model() [2/2]

template<typename PythiaT, typename EventT, typename hepmc_writerT>
void Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::init_user_model ( const std::string  pythiaDocPath,
const std::vector< std::string > &  externalSettings,
std::ostream &  os 
)
inline

Initialize from some external settings, assuming no given SLHAea instance.

Definition at line 197 of file Py8Collider.hpp.

References Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::init_user_model().

199  {
200  init_user_model(pythiaDocPath, externalSettings, nullptr, os);
201  }
void init_user_model(const std::string pythiaDocPath, const std::vector< std::string > &externalSettings, const SLHAea::Coll *slhaea=nullptr, std::ostream &os=std::cout)
Initialize from some external settings.
Here is the call graph for this function:

◆ nextEvent()

template<typename PythiaT, typename EventT, typename hepmc_writerT>
void Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::nextEvent ( EventT &  event) const
inline

Event generation for any Pythia interface to Gambit.

Definition at line 210 of file Py8Collider.hpp.

Referenced by Gambit::ColliderBit::generateEventPy8Collider(), and Gambit::ColliderBit::getPy8Collider().

211  {
212  // Try to make and populate an event
213  bool accepted_event = _pythiaInstance->next();
214  event = _pythiaInstance->event;
215  if (!accepted_event)
216  {
217  throw EventGenerationError();
218  }
219  }
Here is the caller graph for this function:

◆ process_code()

template<typename PythiaT, typename EventT, typename hepmc_writerT>
int Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::process_code ( ) const
inlinevirtual

Report an integer process code for the last generated event.

Implements Gambit::ColliderBit::BaseCollider.

Definition at line 234 of file Py8Collider.hpp.

234 { return _pythiaInstance->info.code(); }

◆ pythia()

template<typename PythiaT, typename EventT, typename hepmc_writerT>
const PythiaT* Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::pythia ( ) const
inline

Get the Pythia instance.

Definition at line 49 of file Py8Collider.hpp.

References Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::_pythiaInstance.

Referenced by Gambit::ColliderBit::generateEventPy8Collider().

49 { return _pythiaInstance; }
Here is the caller graph for this function:

◆ SetupMatchingUserHook()

template<typename PythiaT, typename EventT, typename hepmc_writerT>
bool Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::SetupMatchingUserHook ( )
inline

Definition at line 52 of file Py8Collider.hpp.

References Gambit::ColliderBit::SetHooks< PythiaT, EventT >::SetupHook().

53  {
54  SetHooks<PythiaT, EventT> Hook;
55  Hook.SetupHook(_pythiaInstance);
56  return true;
57  }
Here is the call graph for this function:

◆ xsec_fb() [1/2]

template<typename PythiaT, typename EventT, typename hepmc_writerT>
double Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::xsec_fb ( ) const
inlinevirtual

Report the total or process-specific cross section (in fb or pb).

Implements Gambit::ColliderBit::BaseCollider.

Definition at line 222 of file Py8Collider.hpp.

222 { return _pythiaInstance->info.sigmaGen() * 1e12; }

◆ xsec_fb() [2/2]

template<typename PythiaT, typename EventT, typename hepmc_writerT>
double Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::xsec_fb ( int  process_code) const
inlinevirtual

Implements Gambit::ColliderBit::BaseCollider.

Definition at line 223 of file Py8Collider.hpp.

223 { return _pythiaInstance->info.sigmaGen(process_code) * 1e12; }
int process_code() const
Report an integer process code for the last generated event.

◆ xsec_pb() [1/2]

template<typename PythiaT, typename EventT, typename hepmc_writerT>
double Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::xsec_pb ( ) const
inlinevirtual

Implements Gambit::ColliderBit::BaseCollider.

Definition at line 224 of file Py8Collider.hpp.

224 { return _pythiaInstance->info.sigmaGen() * 1e9; }

◆ xsec_pb() [2/2]

template<typename PythiaT, typename EventT, typename hepmc_writerT>
double Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::xsec_pb ( int  process_code) const
inlinevirtual

Implements Gambit::ColliderBit::BaseCollider.

Definition at line 225 of file Py8Collider.hpp.

225 { return _pythiaInstance->info.sigmaGen(process_code) * 1e9; }
int process_code() const
Report an integer process code for the last generated event.

◆ xsecErr_fb() [1/2]

template<typename PythiaT, typename EventT, typename hepmc_writerT>
double Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::xsecErr_fb ( ) const
inlinevirtual

Report the uncertainty in the total or process-specific cross section (in fb or pb).

Implements Gambit::ColliderBit::BaseCollider.

Definition at line 228 of file Py8Collider.hpp.

228 { return _pythiaInstance->info.sigmaErr() * 1e12; }

◆ xsecErr_fb() [2/2]

template<typename PythiaT, typename EventT, typename hepmc_writerT>
double Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::xsecErr_fb ( int  process_code) const
inlinevirtual

Implements Gambit::ColliderBit::BaseCollider.

Definition at line 229 of file Py8Collider.hpp.

229 { return _pythiaInstance->info.sigmaErr(process_code) * 1e12; }
int process_code() const
Report an integer process code for the last generated event.

◆ xsecErr_pb() [1/2]

template<typename PythiaT, typename EventT, typename hepmc_writerT>
double Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::xsecErr_pb ( ) const
inlinevirtual

Implements Gambit::ColliderBit::BaseCollider.

Definition at line 230 of file Py8Collider.hpp.

230 { return _pythiaInstance->info.sigmaErr() * 1e9; }

◆ xsecErr_pb() [2/2]

template<typename PythiaT, typename EventT, typename hepmc_writerT>
double Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::xsecErr_pb ( int  process_code) const
inlinevirtual

Implements Gambit::ColliderBit::BaseCollider.

Definition at line 231 of file Py8Collider.hpp.

231 { return _pythiaInstance->info.sigmaErr(process_code) * 1e9; }
int process_code() const
Report an integer process code for the last generated event.

Member Data Documentation

◆ _pythiaBase

template<typename PythiaT, typename EventT, typename hepmc_writerT>
PythiaT* Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::_pythiaBase
protected

◆ _pythiaInstance

◆ _pythiaSettings

template<typename PythiaT, typename EventT, typename hepmc_writerT>
std::vector<std::string> Gambit::ColliderBit::Py8Collider< PythiaT, EventT, hepmc_writerT >::_pythiaSettings
protected

Definition at line 43 of file Py8Collider.hpp.


The documentation for this class was generated from the following file: