HepMC3 event record library
ReaderLHEF.cc
Go to the documentation of this file.
1// -*- C++ -*-
2//
3// This file is part of HepMC
4// Copyright (C) 2014-2021 The HepMC collaboration (see AUTHORS for details)
5//
6/**
7 * @file ReaderLHEF.cc
8 * @brief Implementation of \b class ReaderLHEF
9 *
10 */
11#include "HepMC3/ReaderLHEF.h"
12namespace HepMC3
13{
14ReaderLHEF::ReaderLHEF(const std::string& filename)
15{
16 m_reader = std::make_shared<LHEF::Reader>(filename);
17 init();
18}
19ReaderLHEF::ReaderLHEF(std::istream & stream)
20{
21 m_reader = std::make_shared<LHEF::Reader>(stream);
22 init();
23}
24
25ReaderLHEF::ReaderLHEF(std::shared_ptr<std::istream> s_stream)
26 : m_shared_stream(s_stream)
27{
28 m_reader = std::make_shared<LHEF::Reader>(*(m_shared_stream.get()));
29 init();
30}
31
32bool ReaderLHEF::skip(const int n)
33{
34 GenEvent evt;
35 for (int nn = n; nn > 0; --nn)
36 {
37 if (!read_event(evt)) return false;
38 evt.clear();
39 }
40 return !failed();
41}
42
43
45{
46 m_neve = 0;
47 m_failed = false;
48 // Create a HEPRUP attribute and initialize it from the reader.
49 m_hepr = std::make_shared<HEPRUPAttribute>();
50 m_hepr->heprup = m_reader->heprup;
51 // There may be some XML tags in the LHE file which are
52 // non-standard, but we can save them as well.
53 m_hepr->tags = LHEF::XMLTag::findXMLTags(m_reader->headerBlock + m_reader->initComments);
54 // This code is ugly and should be replaced.
55 size_t nweights = 0;
56 for (auto t1: m_hepr->tags) {
57 if (t1->name != "header") continue;
58 for (auto t2: t1->tags) {
59 if (t2->name != "initrwgt") continue;
60 for (auto t3: t2->tags) {
61 if (t3->name != "weightgroup") continue;
62 for (auto t4: t3->tags) if (t4->name == "weight") nweights++;
63 break;
64 }
65 break;
66 }
67 break;
68 }
69 //
70 // Now we want to create a GenRunInfo object for the HepMC file, and
71 // we add the LHEF attribute to that.
72 set_run_info(std::make_shared<GenRunInfo>());
73 run_info()->add_attribute("HEPRUP", m_hepr);
74
75 // This is just a test to make sure we can add other attributes as
76 // well.
77 run_info()->add_attribute("NPRUP",
78 std::make_shared<FloatAttribute>(m_hepr->heprup.NPRUP));
79
80 // We want to be able to convey the different event weights to
81 // HepMC. In particular we need to add the names of the weights to
82 // the GenRunInfo object.
83
84 std::vector<std::string> weightnames;
85 for ( int i = 0, N = m_hepr->heprup.weightinfo.size(); i < N; ++i ) weightnames.push_back(m_hepr->heprup.weightNameHepMC(i));
86 if (nweights == 0) nweights=1;
87 for ( size_t i = weightnames.size(); i < nweights; ++i ) weightnames.push_back(std::to_string(i));
88 run_info()->set_weight_names(weightnames);
89
90 // We also want to convey the information about which generators was
91 // used.
92 for ( int i = 0, N = m_hepr->heprup.generators.size(); i < N; ++i )
93 {
95 tool.name = m_hepr->heprup.generators[i].name;
96 tool.version = m_hepr->heprup.generators[i].version;
97 tool.description = m_hepr->heprup.generators[i].contents;
98 run_info()->tools().push_back(tool);
99 }
100}
101/// @brief Destructor
103
105{
106 if (m_storage.size() > 0)
107 {
108 ev = m_storage.front();
109 m_storage.pop_front();
110 return m_failed;
111 }
112 m_failed = !(m_reader->readEvent());
113 if (m_failed) return m_failed;
114 // To each GenEvent we want to add an attribute corresponding to
115 // the HEPEUP. Also here there may be additional non-standard
116 // information outside the LHEF <event> tags, which we may want to
117 // add.
118 std::shared_ptr<HEPEUPAttribute> hepe = std::make_shared<HEPEUPAttribute>();
119 if ( m_reader->outsideBlock.length() )
120 hepe->tags = LHEF::XMLTag::findXMLTags(m_reader->outsideBlock);
121
122 hepe->hepeup = m_reader->hepeup;
123 std::vector<LHEF::HEPEUP*> input;
124 if (m_reader->hepeup.subevents.size() > 0) input.insert(input.end(), hepe->hepeup.subevents.begin(), hepe->hepeup.subevents.end());
125 else { input.push_back(&m_reader->hepeup);}
126 int first_group_event = m_neve;
127 m_neve++;
128 for (auto ahepeup: input)
129 {
130 GenEvent evt;
131 evt.set_event_number(first_group_event);
132 evt.add_attribute("AlphaQCD", std::make_shared<DoubleAttribute>(ahepeup->AQCDUP));
133 evt.add_attribute("AlphaEM", std::make_shared<DoubleAttribute>(ahepeup->AQEDUP));
134 evt.add_attribute("NUP", std::make_shared<IntAttribute>(ahepeup->NUP));
135 evt.add_attribute("IDPRUP", std::make_shared<LongAttribute>(ahepeup->IDPRUP));
136 // Now add the Particles from the LHE event to HepMC
137 std::vector<GenParticlePtr> particles;
138 std::map< std::pair<int, int>, GenVertexPtr> vertices;
139 for ( int i = 0; i < ahepeup->NUP; ++i )
140 {
141 FourVector mom((ahepeup->PUP)[i][0], (ahepeup->PUP)[i][1], (ahepeup->PUP)[i][2], (ahepeup->PUP)[i][3]);
142 particles.push_back(std::make_shared<GenParticle>(mom, ahepeup->IDUP[i], ahepeup->ISTUP[i]));
143 if ( i < 2 ) continue;
144 std::pair<int, int> vertex_index(ahepeup->MOTHUP[i].first, ahepeup->MOTHUP[i].second);
145 if (vertices.find(vertex_index) == vertices.end()) vertices[vertex_index] = std::make_shared<GenVertex>();
146 vertices[vertex_index]->add_particle_out(particles.back());
147 }
148 for ( auto v: vertices )
149 {
150 std::pair<int, int> vertex_index = v.first;
151 GenVertexPtr vertex = v.second;
152 for (int i = vertex_index.first-1; i < vertex_index.second; ++i)
153 if ( i >= 0 && i < (int)particles.size())
154 vertex->add_particle_in(particles[i]);
155 }
156 std::pair<int, int> vertex_index(0, 0);
157 if (vertices.find(vertex_index) == vertices.end()) vertices[vertex_index] = std::make_shared<GenVertex>();
158 for (size_t i = 0; i < particles.size(); ++i)
159 if (!particles[i]->end_vertex() && !particles[i]->production_vertex())
160 {
161 if ( i < 2 ) vertices[vertex_index]->add_particle_in(particles[i]);
162 else vertices[vertex_index]->add_particle_out(particles[i]);
163 }
164 for ( auto v: vertices ) evt.add_vertex(v.second);
165 if (particles.size() > 1)
166 {
167 particles[0]->set_status(4);
168 particles[1]->set_status(4);
169 evt.set_beam_particles(particles[0], particles[1]);
170 }
171 // And we also want to add the weights.
172
173
174 std::vector<double> wts;
175 for ( int i = 0, N = ahepeup->weights.size(); i < N; ++i )
176 {
177 wts.push_back(ahepeup->weights[i].first);
178 }
179 evt.weights() = wts;
180 m_storage.push_back(evt);
181 }
182 ev = m_storage.front();
183 m_storage.pop_front();
184 return m_failed;
185}
186/// @brief Return status of the stream
188
189/// @brief Close file stream
191} // namespace HepMC3
Definition of class ReaderLHEF.
Generic 4-vector.
Definition: FourVector.h:36
Stores event-related information.
Definition: GenEvent.h:41
void add_vertex(GenVertexPtr v)
Add vertex.
Definition: GenEvent.cc:96
void set_event_number(const int &num)
Set event number.
Definition: GenEvent.h:150
void add_attribute(const std::string &name, const std::shared_ptr< Attribute > &att, const int &id=0)
Definition: GenEvent.cc:806
const std::vector< double > & weights() const
Get event weight values as a vector.
Definition: GenEvent.h:98
void clear()
Remove contents of this event.
Definition: GenEvent.cc:599
void set_beam_particles(GenParticlePtr p1, GenParticlePtr p2)
Set incoming beam particles.
Definition: GenEvent.cc:760
~ReaderLHEF()
Destructor.
Definition: ReaderLHEF.cc:102
void init()
Init helper.
Definition: ReaderLHEF.cc:44
bool failed() override
State.
Definition: ReaderLHEF.cc:187
bool read_event(GenEvent &ev) override
Reading event.
Definition: ReaderLHEF.cc:104
bool skip(const int) override
skip events
Definition: ReaderLHEF.cc:32
std::shared_ptr< HEPRUPAttribute > m_hepr
Holder of attributes.
Definition: ReaderLHEF.h:57
void close() override
Close.
Definition: ReaderLHEF.cc:190
bool m_failed
State of reader.
Definition: ReaderLHEF.h:59
std::deque< GenEvent > m_storage
storage used for subevents.
Definition: ReaderLHEF.h:60
ReaderLHEF(std::istream &)
The ctor to read from stream.
Definition: ReaderLHEF.cc:19
std::shared_ptr< std::istream > m_shared_stream
Holds temporary stream.
Definition: ReaderLHEF.h:55
std::shared_ptr< LHEF::Reader > m_reader
The actual reader.
Definition: ReaderLHEF.h:56
int m_neve
Event counter.
Definition: ReaderLHEF.h:58
std::shared_ptr< GenRunInfo > run_info() const
Get the global GenRunInfo object.
Definition: Reader.h:44
void set_run_info(std::shared_ptr< GenRunInfo > run)
Set the global GenRunInfo object.
Definition: Reader.h:64
HepMC3 main namespace.
Interrnal struct for keeping track of tools.
Definition: GenRunInfo.h:38
std::string description
Other information about how the tool was used in the run.
Definition: GenRunInfo.h:48
std::string version
The version of the tool.
Definition: GenRunInfo.h:44
std::string name
The name of the tool.
Definition: GenRunInfo.h:41
static std::vector< XMLTag * > findXMLTags(std::string str, std::string *leftover=0)
Definition: LHEF.h:198