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