11 #include "Rivet/AnalysisHandler.hh" 12 #include "Rivet/Tools/RivetPaths.hh" 13 #include "Rivet/Config/RivetConfig.hh" 16 #if RIVET_VERSION_CODE - 1 == -1 17 #define RIVET_VERSION_CODE 40000 21 #if RIVET_VERSION_CODE >= 40000 22 #define RIVET_ENABLE_HEPMC_3 23 #define PYTHIA_USING_RIVET_4 27 #ifndef RIVET_ENABLE_HEPMC_3 28 #include "Pythia8Plugins/HepMC2.h" 30 #include "Pythia8Plugins/HepMC3.h" 32 #include "Pythia8/HIInfo.h" 55 : pythia(&pytin), filename(fname), rivet(0), igBeam(false), nDump(-1),
74 void dump(
int period,
string fname =
"") {
81 eventAttributes[name] = val;}
89 #ifndef PYTHIA_USING_RIVET_4 90 rivet =
new Rivet::AnalysisHandler(rname);
91 rivet->setIgnoreBeams(igBeam);
94 rivet->dump(filename, nDump);
96 rivet->dump(dumpFile, nDump);
99 rivet =
new Rivet::AnalysisHandler();
100 rivet->setCheckBeams(!igBeam);
103 rivet->setFinalizePeriod(filename, nDump);
105 rivet->setFinalizePeriod(dumpFile, nDump);
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);
115 rivet->init(converter.
event());
120 for (
const string & path : Rivet::getAnalysisLibPaths() )
121 if ( path ==
"." )
return;
122 Rivet::addAnalysisLibPath(
".");
132 if ( !rivet )
init();
133 #ifdef RIVET_ENABLE_HEPMC_3 134 for (
auto att : eventAttributes)
137 rivet->analyze(converter.
event());
138 eventAttributes.clear();
144 if ( !rivet )
return;
146 rivet->writeData(filename);
160 set<string> analyses;
163 vector<string> preloads;
166 Rivet::AnalysisHandler * rivet;
184 map<const string, double> eventAttributes;
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
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