PYTHIA  8.316
Pythia8Rivet.h
1 // Pythia8Rivet.h is a part of the PYTHIA event generator.
2 // Copyright (C) 2025 Torbjorn Sjostrand.
3 // PYTHIA is licenced under the GNU GPL v2 or later, see COPYING for details.
4 // Please respect the MCnet Guidelines, see GUIDELINES for details.
5 
6 // Author: Christian Bierlich, September 2020.
7 
8 #ifndef PYTHIA8RIVET_H
9 #define PYTHIA8RIVET_H
10 
11 #include "Rivet/AnalysisHandler.hh"
12 #include "Rivet/Tools/RivetPaths.hh"
13 #include "Rivet/Config/RivetConfig.hh"
14 
15 // Check that a proper version code is available (wrong for 4.0.0).
16 #if RIVET_VERSION_CODE - 1 == -1
17 #define RIVET_VERSION_CODE 40000
18 #endif
19 
20 // Deterine which Rivet version is used. If beyond 4, require HepMC3.
21 #if RIVET_VERSION_CODE >= 40000
22 #define RIVET_ENABLE_HEPMC_3
23 #define PYTHIA_USING_RIVET_4
24 #endif
25 
26 // Set the version of HepMC being used (either from Rivet or set above).
27 #ifndef RIVET_ENABLE_HEPMC_3
28 #include "Pythia8Plugins/HepMC2.h"
29 #else
30 #include "Pythia8Plugins/HepMC3.h"
31 #endif
32 #include "Pythia8/HIInfo.h"
33 #include <set>
34 #include <vector>
35 #include <string>
36 
37 namespace Pythia8 {
38 
39 // Simplified interface to the Rivet program. Remember to link with
40 // pythia and -lhepmcinterface -lHepMC -lRivet
41 //
42 // Usage: (1) Create an object giving the pythia object and a filename
43 // as arguments. (2) Repeatedly specify (the name of an) analysis with
44 // the addAnalysis() function, possibly with analysis parameters.
45 // (3) initialize the underlying Rivet object with the init() function.
46 // (4) Analyze an event with the operator() function. (5) Dump the
47 // histograms to a file with the done() function.
48 class Pythia8Rivet {
49 
50 public:
51 
52  // The constructor needs to have the main Pythia object, and the
53  // name of the file where the histograms are dumped.
54  Pythia8Rivet(Pythia & pytin, string fname)
55  : pythia(&pytin), filename(fname), rivet(0), igBeam(false), nDump(-1),
56  dumpFile("") {}
57 
58  // The destructor will write the histogram file if this has not
59  // already been done.
61 
62  // Add the name of an analysis to be performed, with a list
63  // of analysis parameters.
64  void addAnalysis(string ana) {analyses.insert(ana);}
65 
66  // Add a YODA file pre-load pre-filled histogram from.
67  void addPreload(string prel) {preloads.push_back(prel);}
68 
69  // Set "ignore beams" flag.
70  void ignoreBeams(bool flagIn) {igBeam = flagIn;}
71 
72  // Set Rivet dump period and file name for intermittent histograms.
73  // If no dumpfile given, the default output file is used.
74  void dump(int period, string fname = "") {
75  nDump = period;
76  dumpFile = fname;
77  }
78 
79  // Add an attribute to the current event.
80  void addAttribute(const string& name, double val) {
81  eventAttributes[name] = val;}
82 
83  // Add an (optional) run name for Rivet internal use.
84  void addRunName(const string runname) {rname = runname;}
85 
86  // Initialize Rivet. Will do nothing if Rivet was already initialized
87  void init() {
88  if ( rivet ) return;
89 #ifndef PYTHIA_USING_RIVET_4
90  rivet = new Rivet::AnalysisHandler(rname);
91  rivet->setIgnoreBeams(igBeam);
92  if (nDump > 0) {
93  if (dumpFile == "")
94  rivet->dump(filename, nDump);
95  else
96  rivet->dump(dumpFile, nDump);
97  }
98 #else
99  rivet = new Rivet::AnalysisHandler();
100  rivet->setCheckBeams(!igBeam);
101  if (nDump > 0) {
102  if (dumpFile == "")
103  rivet->setFinalizePeriod(filename, nDump);
104  else
105  rivet->setFinalizePeriod(dumpFile, nDump);
106  }
107 #endif
108  initpath();
109  for(int i = 0, N = preloads.size(); i < N; ++i)
110  rivet->readData(preloads[i]);
111  for (set<string>::iterator it = analyses.begin();
112  it != analyses.end(); ++it) {
113  rivet->addAnalysis(*it);
114  }
115  rivet->init(converter.event());
116  }
117 
118  // Helper function to set analysis path.
119  void initpath() const {
120  for ( const string & path : Rivet::getAnalysisLibPaths() )
121  if ( path == "." ) return;
122  Rivet::addAnalysisLibPath(".");
123  }
124 
125  // Analyze the default Pythia event. Will do nothing if Rivet has
126  // not been intialized.
127  void operator()() {this->operator()(pythia->event);}
128 
129  // Analyze the given event.
130  void operator()(Event & event) {
131  converter.fillNextEvent(*pythia);
132  if ( !rivet ) init();
133 #ifdef RIVET_ENABLE_HEPMC_3
134  for (auto att : eventAttributes)
135  converter.addAttribute(att.first, att.second);
136 #endif
137  rivet->analyze(converter.event());
138  eventAttributes.clear();
139  }
140 
141  // Writes histograms to file and deletes the Rivet object. Does
142  // nothing if Rivet was not initialized.
143  void done() {
144  if ( !rivet ) return;
145  rivet->finalize();
146  rivet->writeData(filename);
147  delete rivet;
148  rivet = 0;
149  }
150 
151 private:
152 
153  // The main pythia object.
154  Pythia * pythia;
155 
156  // The name of the file where the histograms are dumped.
157  string filename;
158 
159  // Analyses with optional analysis parameters.
160  set<string> analyses;
161 
162  // The names of YODA files to preload.
163  vector<string> preloads;
164 
165  // The Rivet object.
166  Rivet::AnalysisHandler * rivet;
167 
168  // The HepMC converter.
169  Pythia8ToHepMC converter;
170 
171  // The Rivet run name.
172  string rname;
173 
174  // Ignore beams flag.
175  bool igBeam;
176 
177  // Dump period.
178  int nDump;
179 
180  // Optional alternative dumpfile name.
181  string dumpFile;
182 
183  // Additional event attributes of double type to send to Rivet.
184  map<const string, double> eventAttributes;
185 
186 };
187 
188 }
189 
190 #endif
GenEvent & event()
Get a reference to the current GenEvent.
Definition: HepMC2.h:465
The Event class holds all info on the generated event.
Definition: Event.h:408
void init()
Initialize Rivet. Will do nothing if Rivet was already initialized.
Definition: Pythia8Rivet.h:87
Event event
The event record for the complete event history.
Definition: Pythia.h:375
void addAttribute(const string &name, T &attribute)
Definition: HepMC3.h:421
void addAnalysis(string ana)
Definition: Pythia8Rivet.h:64
~Pythia8Rivet()
Definition: Pythia8Rivet.h:60
void dump(int period, string fname="")
Definition: Pythia8Rivet.h:74
void done()
Definition: Pythia8Rivet.h:143
Definition: HepMC2.h:409
void addAttribute(const string &name, double val)
Add an attribute to the current event.
Definition: Pythia8Rivet.h:80
void operator()(Event &event)
Analyze the given event.
Definition: Pythia8Rivet.h:130
void operator()()
Definition: Pythia8Rivet.h:127
void ignoreBeams(bool flagIn)
Set "ignore beams" flag.
Definition: Pythia8Rivet.h:70
Pythia8Rivet(Pythia &pytin, string fname)
Definition: Pythia8Rivet.h:54
void addRunName(const string runname)
Add an (optional) run name for Rivet internal use.
Definition: Pythia8Rivet.h:84
Header for classes to set beam momentum and interaction vertex spread.
Definition: Analysis.h:20
The Pythia class contains the top-level routines to generate an event.
Definition: Pythia.h:72
bool fillNextEvent(Pythia &pythia)
Definition: HepMC2.h:445
void addPreload(string prel)
Add a YODA file pre-load pre-filled histogram from.
Definition: Pythia8Rivet.h:67
Definition: Pythia8Rivet.h:48
void initpath() const
Helper function to set analysis path.
Definition: Pythia8Rivet.h:119