HepMC3 event record library
LHEF_example_cat.cc
1/**
2 * @example LHEF_example_cat.cc
3 * @brief Basic example of use of LHEF for reading and writing LHE files
4 */
8#include "HepMC3/GenEvent.h"
10#include "HepMC3/GenVertex.h"
12#include <iomanip>
13
14using namespace HepMC3;
15
16int main(int /*argc*/, char ** /*argv*/) {
17
18 // Create Reader to read the example LHE file.
19 LHEF::Reader reader("LHEF_example.lhe");
20
21 // Create a HEPRUP attribute and initialize it from the reader.
22 std::shared_ptr<HEPRUPAttribute> hepr = std::make_shared<HEPRUPAttribute>();
23 hepr->heprup = reader.heprup;
24
25 // There may be some XML tags in the LHE file which are
26 // non-standard, but we can save them as well.
27 hepr->tags = LHEF::XMLTag::findXMLTags(reader.headerBlock + reader.initComments);
28
29 // Nowwe want to create a GenRunInfo object for the HepMC file, and
30 // we add the LHEF attribute to that.
31 std::shared_ptr<GenRunInfo> runinfo = std::make_shared<GenRunInfo>();
32 runinfo->add_attribute("HEPRUP", hepr);
33
34 // This is just a test to make sure we can add other attributes as
35 // well.
36 runinfo->add_attribute("NPRUP",
37 std::make_shared<FloatAttribute>(hepr->heprup.NPRUP));
38
39 // We want to be able to convey the different event weights to
40 // HepMC. In particular we need to add the names of the weights to
41 // the GenRunInfo object.
42 std::vector<std::string> weightnames;
43 weightnames.emplace_back("0"); // The first weight is always the
44 // default weight with name "0".
45 for ( int i = 0, N = hepr->heprup.weightinfo.size(); i < N; ++i ){
46 weightnames.push_back(hepr->heprup.weightNameHepMC(i));}
47 runinfo->set_weight_names(weightnames);
48
49 // We also want to convey the information about which generators was
50 // used to HepMC.
51 for ( int i = 0, N = hepr->heprup.generators.size(); i < N; ++i ) {
53 tool.name = hepr->heprup.generators[i].name;
54 tool.version = hepr->heprup.generators[i].version;
55 tool.description = hepr->heprup.generators[i].contents;
56 runinfo->tools().push_back(tool);
57 }
58
59 // Now we want to start reading events from the LHE file and
60 // translate them to HepMC.
61 WriterAscii output("LHEF_example.hepmc3", runinfo);
62 int neve = 0;
63 while ( reader.readEvent() ) {
64 ++neve;
65
66 // To each GenEvent we want to add an attribute corresponding to
67 // the HEPEUP. Also here there may be additional non-standard
68 // information outside the LHEF <event> tags, which we may want to
69 // add.
70 std::shared_ptr<HEPEUPAttribute> hepe = std::make_shared<HEPEUPAttribute>();
71 if ( reader.outsideBlock.length() ){
72 hepe->tags = LHEF:: XMLTag::findXMLTags(reader.outsideBlock);
73 }
74 hepe->hepeup = reader.hepeup;
75 GenEvent ev(runinfo, Units::GEV, Units::MM);
76 ev.set_event_number(neve);
77
78 // This is just a text to check that we can add additional
79 // attributes to each event.
80 ev.add_attribute("HEPEUP", hepe);
81 ev.add_attribute("AlphaQCD",
82 std:: make_shared<DoubleAttribute>(hepe->hepeup.AQCDUP));
83 ev.add_attribute("AlphaEM",
84 std::make_shared<DoubleAttribute>(hepe->hepeup.AQEDUP));
85 ev.add_attribute("NUP",
86 std::make_shared<IntAttribute>(hepe->hepeup.NUP));
87 ev.add_attribute("IDPRUP",
88 std::make_shared<LongAttribute>(hepe->hepeup.IDPRUP));
89
90 // Now add the Particles from the LHE event to HepMC
91 std::vector<GenParticlePtr> particles;
92 std::map< std::pair<int,int>, GenVertexPtr> vertices;
93 for ( int i = 0; i < hepe->hepeup.NUP; ++i )
94 {
95 particles.push_back(std::make_shared<GenParticle>(hepe->momentum(i),hepe->hepeup.IDUP[i],hepe->hepeup.ISTUP[i]));
96 if (i<2) continue;
97 std::pair<int,int> vertex_index(hepe->hepeup.MOTHUP[i].first,hepe->hepeup.MOTHUP[i].second);
98 if (vertices.find(vertex_index)==vertices.end()) {vertices[vertex_index]=std::make_shared<GenVertex>();}
99 vertices[vertex_index]->add_particle_out(particles.back());
100 }
101 for ( auto& v: vertices )
102 {
103 std::pair<int,int> vertex_index=v.first;
104 GenVertexPtr vertex=v.second;
105 for (int i=vertex_index.first-1; i<vertex_index.second; i++) {if (i >= 0 && i < static_cast<int>(particles.size())) vertex->add_particle_in(particles[i]);}
106 }
107 for ( auto& v: vertices ) {ev.add_vertex(v.second);}
108
109 // And we also want to add the weights.
110 std::vector<double> wts;
111 for ( int i = 0, N = hepe->hepeup.weights.size(); i < N; ++i ){
112 wts.push_back(hepe->hepeup.weights[i].first);
113 }
114 ev.weights() = wts;
115
116 // Let's see if we can associate p1 and p2.
117 ev.add_attribute("OtherIncoming",
118 std::make_shared<AssociatedParticle>(particles[1]), particles[0]->id());
119
120
121 // And then we are done and can write out the GenEvent.
122 output.write_event(ev);
123
124 }
125
126 output.close();
127
128 // Now we wnat to make sure we can read in the HepMC file and
129 // recreate the same info. To check this we try to reconstruct the
130 // LHC file we read in above.
131 ReaderAscii input("LHEF_example.hepmc3");
132 LHEF::Writer writer("LHEF_example_out.lhe");
133 hepr = std::shared_ptr<HEPRUPAttribute>();
134
135 // The loop over all events in the HepMC file.
136 while ( true ) {
137
138 // Read in the next event.
139 GenEvent ev(Units::GEV, Units::MM);
140 if ( !input.read_event(ev) || ev.event_number() == 0 ) break;
141
142 // Check that the first incoming particle still refers to the second.
143 std::shared_ptr<AssociatedParticle> assoc =
144 ev.attribute<AssociatedParticle>("OtherIncoming", 1);
145 if ( !assoc || !assoc->associated() ||
146 assoc->associated() != ev.particles()[1] ) return 3;
147
148 // Make sure the weight names are the same.
149 if ( input.run_info()->weight_names() != weightnames ) return 2;
150
151 // For the first event we also go in and reconstruct the HEPRUP
152 // information, and write it out to the new LHE file.
153 if ( !hepr ) {
154 hepr = ev.attribute<HEPRUPAttribute>("HEPRUP");
155
156 // Here we also keep track of the additional non-standard info
157 // we found in the original LHE file.
158 for ( int i = 0, N = hepr->tags.size(); i < N; ++i ){
159 if ( hepr->tags[i]->name != "init" ){
160 hepr->tags[i]->print(writer.headerBlock());
161 }
162}
163 // This is just a test that we can access other attributes
164 // included in the GenRunInfo.
165 hepr->heprup.NPRUP =
166 int(input.run_info()->attribute<FloatAttribute>("NPRUP")->value());
167
168 // Then we write out the HEPRUP object.
169 writer.heprup = hepr->heprup;
170 if ( writer.heprup.eventfiles.size() >= 2 ) {
171 writer.heprup.eventfiles[0].filename = "LHEF_example_1_out.plhe";
172 writer.heprup.eventfiles[1].filename = "LHEF_example_2_out.plhe";
173 }
174 writer.init();
175
176 }
177
178 // Now we can access the HEPEUP attribute of the current event.
179 std::shared_ptr<HEPEUPAttribute> hepe = ev.attribute<HEPEUPAttribute>("HEPEUP");
180
181 // Again, there may be addisional non-standard information we want
182 // to keep.
183 for ( int i = 0, N = hepe->tags.size(); i < N; ++i ){
184 if ( hepe->tags[i]->name != "event" &&
185 hepe->tags[i]->name != "eventgroup" ){
186 hepe->tags[i]->print(writer.eventComments());
187}
188}
189 // This is just a test that we can access other attributes
190 // included in the GenRunInfo.
191 hepe->hepeup.AQCDUP =
192 ev.attribute<DoubleAttribute>("AlphaQCD")->value();
193 hepe->hepeup.AQEDUP =
194 ev.attribute<DoubleAttribute>("AlphaEM")->value();
195 hepe->hepeup.NUP =
196 ev.attribute<IntAttribute>("NUP")->value();
197 hepe->hepeup.IDPRUP =
198 ev.attribute<LongAttribute>("IDPRUP")->value();
199
200 // And then we can write out the HEPEUP object.
201 writer.hepeup = hepe->hepeup;
202 writer.hepeup.heprup = &writer.heprup;
203 writer.writeEvent();
204
205 }
206
207}
Definition of class AssociatedParticle,.
Definition of class GenEvent.
Definition of class GenParticle.
Definition of class GenVertex.
Definition of class HEPRUPAttribute and class HEPEUAttribute.
Definition of class ReaderAscii.
Definition of class WriterAscii.
Attribute class allowing eg. a GenParticle to refer to another GenParticle.
Attribute that holds a real number as a double.
Definition Attribute.h:245
Attribute that holds a real number as a float.
Definition Attribute.h:292
float value() const
get the value associated to this Attribute.
Definition Attribute.h:318
Stores event-related information.
Definition GenEvent.h:47
Class for storing data for LHEF run information.
Class for storing data for LHEF run information.
Attribute that holds an Integer implemented as an int.
Definition Attribute.h:157
Attribute that holds an Integer implemented as a long int.
Definition Attribute.h:200
GenEvent I/O parsing for structured text files.
Definition ReaderAscii.h:31
GenEvent I/O serialization for structured text files.
Definition WriterAscii.h:25
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=nullptr)
Definition LHEF.h:200