HepMC3 event record library
GenEvent.h
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // This file is part of HepMC
4 // Copyright (C) 2014-2019 The HepMC collaboration (see AUTHORS for details)
5 //
6 ///
7 /// @file GenEvent.h
8 /// @brief Definition of \b class GenEvent
9 ///
10 #ifndef HEPMC3_GENEVENT_H
11 #define HEPMC3_GENEVENT_H
12 
13 #include "HepMC3/Units.h"
14 #include "HepMC3/GenParticle_fwd.h"
15 #include "HepMC3/GenVertex_fwd.h"
16 #include "HepMC3/GenPdfInfo_fwd.h"
17 #include "HepMC3/GenHeavyIon_fwd.h"
18 #include "HepMC3/GenCrossSection_fwd.h"
19 
20 #if !defined(__CINT__)
21 #include "HepMC3/Errors.h"
22 #include "HepMC3/GenHeavyIon.h"
23 #include "HepMC3/GenPdfInfo.h"
24 #include "HepMC3/GenCrossSection.h"
25 #include "HepMC3/GenRunInfo.h"
26 #include <mutex>
27 #endif // __CINT__
28 
29 #ifdef HEPMC3_ROOTIO
30 class TBuffer;
31 #endif
32 
33 
34 namespace HepMC3 {
35 
36 struct GenEventData;
37 
38 /// @brief Stores event-related information
39 ///
40 /// Manages event-related information.
41 /// Contains lists of GenParticle and GenVertex objects
42 class GenEvent {
43 
44 public:
45 
46  /// @brief Event constructor without a run
48  Units::LengthUnit length_unit = Units::MM);
49 
50 #if !defined(__CINT__)
51 
52  /// @brief Constructor with associated run
53  GenEvent(shared_ptr<GenRunInfo> run,
55  Units::LengthUnit length_unit = Units::MM);
56 
57  /// @brief Copy constructor
58  GenEvent(const GenEvent&);
59 
60  /// @brief Destructor
61  ~GenEvent();
62 
63  /// @brief Assignment operator
64  GenEvent& operator=(const GenEvent&);
65 
66  /// @name Particle and vertex access
67  //@{
68 
69  /// @brief Get list of particles (const)
70  const std::vector<ConstGenParticlePtr>& particles() const;
71  /// @brief Get list of vertices (const)
72  const std::vector<ConstGenVertexPtr>& vertices() const;
73 
74 
75  /// @brief Get/set list of particles (non-const)
76  const std::vector<GenParticlePtr>& particles() { return m_particles; }
77  /// @brief Get/set list of vertices (non-const)
78  const std::vector<GenVertexPtr>& vertices() { return m_vertices; }
79 
80  //@}
81 
82 
83  /// @name Event weights
84  //@{
85 
86  /// Get event weight values as a vector
87  const std::vector<double>& weights() const { return m_weights; }
88  /// Get event weights as a vector (non-const)
89  std::vector<double>& weights() { return m_weights; }
90  /// Get event weight accessed by index (or the canonical/first one if there is no argument)
91  /// @note It's the user's responsibility to ensure that the given index exists!
92  double weight(const size_t& index=0) const { return weights().at(index); }
93  /// Get event weight accessed by weight name
94  /// @note Requires there to be an attached GenRunInfo, otherwise will throw an exception
95  /// @note It's the user's responsibility to ensure that the given name exists!
96  double weight(const std::string& name) const {
97  if (!run_info()) throw WeightError("GenEvent::weight(str): named access to event weights requires the event to have a GenRunInfo");
98  return weight(run_info()->weight_index(name));
99  }
100  /// Get event weight accessed by weight name
101  /// @note Requires there to be an attached GenRunInfo, otherwise will throw an exception
102  /// @note It's the user's responsibility to ensure that the given name exists!
103  double& weight(const std::string& name) {
104  if (!run_info()) throw WeightError("GenEvent::weight(str): named access to event weights requires the event to have a GenRunInfo");
105  int pos=run_info()->weight_index(name);
106  if (pos<0) throw WeightError("GenEvent::weight(str): no weight with given name in this run");
107  return m_weights[pos];
108  }
109  /// Get event weight names, if there are some
110  /// @note Requires there to be an attached GenRunInfo with registered weight names, otherwise will throw an exception
111  const std::vector<std::string>& weight_names(const std::string& /*name*/) const {
112  if (!run_info()) throw WeightError("GenEvent::weight_names(): access to event weight names requires the event to have a GenRunInfo");
113  const std::vector<std::string>& weightnames = run_info()->weight_names();
114  if (weightnames.empty()) throw WeightError("GenEvent::weight_names(): no event weight names are registered for this run");
115  return weightnames;
116  }
117 
118  //@}
119 
120 
121  /// @name Auxiliary info and event metadata
122  //@{
123 
124  /// @brief Get a pointer to the the GenRunInfo object.
125  shared_ptr<GenRunInfo> run_info() const {
126  return m_run_info;
127  }
128  /// @brief Set the GenRunInfo object by smart pointer.
129  void set_run_info(shared_ptr<GenRunInfo> run) {
130  m_run_info = run;
131  if ( run && !run->weight_names().empty() )
132  m_weights.resize(run->weight_names().size(), 1.0);
133  }
134 
135  /// @brief Get event number
136  int event_number() const { return m_event_number; }
137  /// @brief Set event number
138  void set_event_number(const int& num) { m_event_number = num; }
139 
140  /// @brief Get momentum unit
142  /// @brief Get length unit
143  const Units::LengthUnit& length_unit() const { return m_length_unit; }
144  /// @brief Change event units
145  /// Converts event from current units to new ones
146  void set_units( Units::MomentumUnit new_momentum_unit, Units::LengthUnit new_length_unit);
147 
148  /// @brief Get heavy ion generator additional information
149  GenHeavyIonPtr heavy_ion() { return attribute<GenHeavyIon>("GenHeavyIon"); }
150  /// @brief Get heavy ion generator additional information (const version)
151  ConstGenHeavyIonPtr heavy_ion() const { return attribute<GenHeavyIon>("GenHeavyIon"); }
152  /// @brief Set heavy ion generator additional information
153  void set_heavy_ion(GenHeavyIonPtr hi) { add_attribute("GenHeavyIon",hi); }
154 
155  /// @brief Get PDF information
156  GenPdfInfoPtr pdf_info() { return attribute<GenPdfInfo>("GenPdfInfo"); }
157  /// @brief Get PDF information (const version)
158  ConstGenPdfInfoPtr pdf_info() const { return attribute<GenPdfInfo>("GenPdfInfo"); }
159  /// @brief Set PDF information
160  void set_pdf_info(GenPdfInfoPtr pi) { add_attribute("GenPdfInfo",pi); }
161 
162  /// @brief Get cross-section information
163  GenCrossSectionPtr cross_section() { return attribute<GenCrossSection>("GenCrossSection"); }
164  /// @brief Get cross-section information (const version)
165  ConstGenCrossSectionPtr cross_section() const { return attribute<GenCrossSection>("GenCrossSection"); }
166  /// @brief Set cross-section information
167  void set_cross_section(GenCrossSectionPtr cs) { add_attribute("GenCrossSection",cs); }
168 
169  //@}
170 
171 
172  /// @name Event position
173  //@{
174 
175  /// Vertex representing the overall event position
176  const FourVector& event_pos() const;
177 
178  /// @brief Vector of beam particles
179  std::vector<ConstGenParticlePtr> beams() const;
180 
181  /// @brief Vector of beam particles
182  const std::vector<GenParticlePtr> & beams();
183 
184  /// @brief Shift position of all vertices in the event by @a delta
185  void shift_position_by( const FourVector & delta );
186 
187  /// @brief Shift position of all vertices in the event to @a op
188  void shift_position_to( const FourVector & newpos ) {
189  const FourVector delta = newpos - event_pos();
190  shift_position_by(delta);
191  }
192 
193  /// @brief Boost event using x,y,z components of @a v as velocities
194  bool boost( const FourVector v );
195  /// @brief Rotate event using x,y,z components of @a v as rotation angles
196  bool rotate( const FourVector v );
197  /// @brief Change sign of @a axis
198  bool reflect(const int axis);
199 
200  //@}
201 
202 
203  /// @name Additional attributes
204  //@{
205  /// @brief Add event attribute to event
206  ///
207  /// This will overwrite existing attribute if an attribute
208  /// with the same name is present
209  void add_attribute(const string &name, const shared_ptr<Attribute> &att, const int& id = 0) {
210  std::lock_guard<std::recursive_mutex> lock(m_lock_attributes);
211  if ( att ) {
212  m_attributes[name][id] = att;
213  att->m_event = this;
214  if ( id > 0 && id <= int(particles().size()) )
215  att->m_particle = particles()[id - 1];
216  if ( id < 0 && -id <= int(vertices().size()) )
217  att->m_vertex = vertices()[-id - 1];
218  }
219  }
220 
221  /// @brief Remove attribute
222  void remove_attribute(const string &name, const int& id = 0);
223 
224  /// @brief Get attribute of type T
225  template<class T>
226  shared_ptr<T> attribute(const string &name, const int& id = 0) const;
227 
228  /// @brief Get attribute of any type as string
229  string attribute_as_string(const string &name, const int& id = 0) const;
230 
231  /// @brief Get list of attribute names
232  std::vector<string> attribute_names( const int& id = 0) const;
233 
234  /// @brief Get a copy of the list of attributes
235  /// @note To avoid thread issues, this is returns a copy. Better solution may be needed.
236  std::map< string, std::map<int, shared_ptr<Attribute> > > attributes() const {
237  std::lock_guard<std::recursive_mutex> lock(m_lock_attributes);
238  return m_attributes;
239  }
240 
241  //@}
242 
243 
244  /// @name Particle and vertex modification
245  //@{
246 
247  /// @brief Add particle
248  void add_particle( GenParticlePtr p );
249 
250  /// @brief Add vertex
251  void add_vertex( GenVertexPtr v );
252 
253  /// @brief Remove particle from the event
254  ///
255  /// This function will remove whole sub-tree starting from this particle
256  /// if it is the only incoming particle of this vertex.
257  /// It will also production vertex of this particle if this vertex
258  /// has no more outgoing particles
259  void remove_particle( GenParticlePtr v );
260 
261  /// @brief Remove a set of particles
262  ///
263  /// This function follows rules of GenEvent::remove_particle to remove
264  /// a list of particles from the event.
265  void remove_particles( std::vector<GenParticlePtr> v );
266 
267  /// @brief Remove vertex from the event
268  ///
269  /// This will remove all sub-trees of all outgoing particles of this vertex
270  void remove_vertex( GenVertexPtr v );
271 
272  /// @brief Add whole tree in topological order
273  ///
274  /// This function will find the beam particles (particles
275  /// that have no production vertices or their production vertices
276  /// have no particles) and will add the whole decay tree starting from
277  /// these particles.
278  ///
279  /// @note Any particles on this list that do not belong to the tree
280  /// will be ignored.
281  void add_tree( const std::vector<GenParticlePtr> &particles );
282 
283  /// @brief Reserve memory for particles and vertices
284  ///
285  /// Helps optimize event creation when size of the event is known beforehand
286  void reserve(const size_t& particles, const size_t& vertices = 0);
287 
288  /// @brief Remove contents of this event
289  void clear();
290 
291  //@}
292 
293  /// @name Deprecated functionality
294  //@{
295 
296  /// @brief Add particle by raw pointer
297  /// @deprecated Use GenEvent::add_particle( const GenParticlePtr& ) instead
298  void add_particle( GenParticle *p );
299 
300  /// @brief Add vertex by raw pointer
301  /// @deprecated Use GenEvent::add_vertex( const GenVertexPtr& ) instead
302  void add_vertex ( GenVertex *v );
303 
304 
305  /// @brief Set incoming beam particles
306  /// @deprecated Backward compatibility
307  void set_beam_particles(GenParticlePtr p1, GenParticlePtr p2);
308 
309 
310  /// @brief Add particle to root vertex
311 
312  void add_beam_particle(GenParticlePtr p1);
313 
314 
315  //@}
316 
317 #endif // __CINT__
318 
319 
320  /// @name Methods to fill GenEventData and to read it back
321  //@{
322 
323  /// @brief Fill GenEventData object
324  void write_data(GenEventData &data) const;
325 
326  /// @brief Fill GenEvent based on GenEventData
327  void read_data(const GenEventData &data);
328 
329 #ifdef HEPMC3_ROOTIO
330  /// @brief ROOT I/O streamer
331  void Streamer(TBuffer &b);
332  //@}
333 #endif
334 
335 private:
336 
337  /// @name Fields
338  //@{
339 
340 #if !defined(__CINT__)
341 
342  /// List of particles
343  std::vector<GenParticlePtr> m_particles;
344  /// List of vertices
345  std::vector<GenVertexPtr> m_vertices;
346 
347  /// Event number
349 
350  /// Event weights
351  std::vector<double> m_weights;
352 
353  /// Momentum unit
355  /// Length unit
357 
358  /// The root vertex is stored outside the normal vertices list to block user access to it
359  GenVertexPtr m_rootvertex;
360 
361  /// Global run information.
362  shared_ptr<GenRunInfo> m_run_info;
363 
364  /// @brief Map of event, particle and vertex attributes
365  ///
366  /// Keys are name and ID (0 = event, <0 = vertex, >0 = particle)
367  mutable std::map< string, std::map<int, shared_ptr<Attribute> > > m_attributes;
368 
369  /// @brief Attribute map key type
370  typedef std::map< string, std::map<int, shared_ptr<Attribute> > >::value_type att_key_t;
371 
372  /// @brief Attribute map value type
373  typedef std::map<int, shared_ptr<Attribute> >::value_type att_val_t;
374 
375  /// @brief Mutex lock for the m_attibutes map.
376  mutable std::recursive_mutex m_lock_attributes;
377 #endif // __CINT__
378 
379  //@}
380 
381 };
382 
383 #if !defined(__CINT__)
384 //
385 // Template methods
386 //
387 template<class T>
388 shared_ptr<T> GenEvent::attribute(const std::string &name, const int& id) const {
389  std::lock_guard<std::recursive_mutex> lock(m_lock_attributes);
390  std::map< string, std::map<int, shared_ptr<Attribute> > >::iterator i1 =
391  m_attributes.find(name);
392  if( i1 == m_attributes.end() ) {
393  if ( id == 0 && run_info() ) {
394  return run_info()->attribute<T>(name);
395  }
396  return shared_ptr<T>();
397  }
398 
399  std::map<int, shared_ptr<Attribute> >::iterator i2 = i1->second.find(id);
400  if (i2 == i1->second.end() ) return shared_ptr<T>();
401 
402  if (!i2->second->is_parsed() ) {
403 
404  shared_ptr<T> att = make_shared<T>();
405  att->m_event = this;
406 
407  if ( id > 0 && id <= int(particles().size()) )
408  att->m_particle = m_particles[id - 1];
409  if ( id < 0 && -id <= int(vertices().size()) )
410  att->m_vertex = m_vertices[-id - 1];
411  if ( att->from_string(i2->second->unparsed_string()) &&
412  att->init() ) {
413  // update map with new pointer
414  i2->second = att;
415  return att;
416  } else {
417  return shared_ptr<T>();
418  }
419  }
420  else return dynamic_pointer_cast<T>(i2->second);
421 }
422 #endif // __CINT__
423 
424 } // namespace HepMC3
425 #endif
HepMC3::GenEvent::heavy_ion
GenHeavyIonPtr heavy_ion()
Get heavy ion generator additional information.
Definition: GenEvent.h:149
HepMC3::FourVector
Generic 4-vector.
Definition: FourVector.h:35
HepMC3::GenEvent::event_number
int event_number() const
Get event number.
Definition: GenEvent.h:136
HepMC3::GenEvent::operator=
GenEvent & operator=(const GenEvent &)
Assignment operator.
Definition: GenEvent.cc:84
HepMC3::WeightError
Exception related to weight lookups, setting, and index consistency.
Definition: Errors.h:56
HepMC3::GenEvent::set_pdf_info
void set_pdf_info(GenPdfInfoPtr pi)
Set PDF information.
Definition: GenEvent.h:160
HepMC3::GenEvent::attributes
std::map< string, std::map< int, shared_ptr< Attribute > > > attributes() const
Get a copy of the list of attributes.
Definition: GenEvent.h:236
HepMC3::GenEvent
Stores event-related information.
Definition: GenEvent.h:42
HepMC3::GenEvent::add_vertex
void add_vertex(GenVertexPtr v)
Add vertex.
Definition: GenEvent.cc:98
HepMC3::GenEvent::reflect
bool reflect(const int axis)
Change sign of axis.
Definition: GenEvent.cc:527
HepMC3::GenEvent::particles
const std::vector< GenParticlePtr > & particles()
Get/set list of particles (non-const)
Definition: GenEvent.h:76
HepMC3::GenEvent::boost
bool boost(const FourVector v)
Boost event using x,y,z components of v as velocities.
Definition: GenEvent.cc:559
HepMC3::GenEvent::cross_section
ConstGenCrossSectionPtr cross_section() const
Get cross-section information (const version)
Definition: GenEvent.h:165
HepMC3
HepMC3 main namespace.
Definition: ReaderGZ.h:28
HepMC3::GenEvent::set_beam_particles
void set_beam_particles(GenParticlePtr p1, GenParticlePtr p2)
Set incoming beam particles.
Definition: GenEvent.cc:776
HepMC3::GenEventData
Stores serializable event information.
Definition: GenEventData.h:26
HepMC3::GenEvent::pdf_info
ConstGenPdfInfoPtr pdf_info() const
Get PDF information (const version)
Definition: GenEvent.h:158
HepMC3::Units::MomentumUnit
MomentumUnit
Momentum units.
Definition: Units.h:29
HepMC3::GenEvent::beams
std::vector< ConstGenParticlePtr > beams() const
Vector of beam particles.
Definition: GenEvent.cc:421
HepMC3::GenEvent::cross_section
GenCrossSectionPtr cross_section()
Get cross-section information.
Definition: GenEvent.h:163
HepMC3::GenEvent::clear
void clear()
Remove contents of this event.
Definition: GenEvent.cc:610
HepMC3::GenEvent::momentum_unit
const Units::MomentumUnit & momentum_unit() const
Get momentum unit.
Definition: GenEvent.h:141
HepMC3::GenEvent::pdf_info
GenPdfInfoPtr pdf_info()
Get PDF information.
Definition: GenEvent.h:156
HepMC3::GenEvent::reserve
void reserve(const size_t &particles, const size_t &vertices=0)
Reserve memory for particles and vertices.
Definition: GenEvent.cc:390
HepMC3::GenEvent::vertices
const std::vector< GenVertexPtr > & vertices()
Get/set list of vertices (non-const)
Definition: GenEvent.h:78
HepMC3::GenEvent::shift_position_by
void shift_position_by(const FourVector &delta)
Shift position of all vertices in the event by delta.
Definition: GenEvent.cc:429
HepMC3::GenEvent::set_run_info
void set_run_info(shared_ptr< GenRunInfo > run)
Set the GenRunInfo object by smart pointer.
Definition: GenEvent.h:129
HepMC3::GenEvent::vertices
const std::vector< ConstGenVertexPtr > & vertices() const
Get list of vertices (const)
Definition: GenEvent.cc:44
GenHeavyIon.h
Definition of attribute class GenHeavyIon.
HepMC3::GenEvent::heavy_ion
ConstGenHeavyIonPtr heavy_ion() const
Get heavy ion generator additional information (const version)
Definition: GenEvent.h:151
HepMC3::GenEvent::write_data
void write_data(GenEventData &data) const
Fill GenEventData object.
Definition: GenEvent.cc:647
GenPdfInfo.h
Definition of event attribute class GenPdfInfo.
HepMC3::GenEvent::m_event_number
int m_event_number
Event number.
Definition: GenEvent.h:348
HepMC3::GenEvent::weights
std::vector< double > & weights()
Get event weights as a vector (non-const)
Definition: GenEvent.h:89
HepMC3::GenEvent::m_momentum_unit
Units::MomentumUnit m_momentum_unit
Momentum unit.
Definition: GenEvent.h:354
HepMC3::GenEvent::set_event_number
void set_event_number(const int &num)
Set event number.
Definition: GenEvent.h:138
HepMC3::GenEvent::weights
const std::vector< double > & weights() const
Get event weight values as a vector.
Definition: GenEvent.h:87
HepMC3::GenEvent::event_pos
const FourVector & event_pos() const
Vertex representing the overall event position.
Definition: GenEvent.cc:417
HepMC3::GenEvent::shift_position_to
void shift_position_to(const FourVector &newpos)
Shift position of all vertices in the event to op.
Definition: GenEvent.h:188
HepMC3::GenEvent::m_run_info
shared_ptr< GenRunInfo > m_run_info
Global run information.
Definition: GenEvent.h:362
HepMC3::GenEvent::run_info
shared_ptr< GenRunInfo > run_info() const
Get a pointer to the the GenRunInfo object.
Definition: GenEvent.h:125
HepMC3::GenEvent::add_beam_particle
void add_beam_particle(GenParticlePtr p1)
Add particle to root vertex.
Definition: GenEvent.cc:781
HepMC3::GenEvent::m_lock_attributes
std::recursive_mutex m_lock_attributes
Mutex lock for the m_attibutes map.
Definition: GenEvent.h:376
HepMC3::GenEvent::rotate
bool rotate(const FourVector v)
Rotate event using x,y,z components of v as rotation angles.
Definition: GenEvent.cc:439
HepMC3::GenEvent::weight
double & weight(const std::string &name)
Definition: GenEvent.h:103
HepMC3::GenEvent::m_weights
std::vector< double > m_weights
Event weights.
Definition: GenEvent.h:351
HepMC3::GenEvent::add_attribute
void add_attribute(const string &name, const shared_ptr< Attribute > &att, const int &id=0)
Add event attribute to event.
Definition: GenEvent.h:209
HepMC3::GenEvent::remove_attribute
void remove_attribute(const string &name, const int &id=0)
Remove attribute.
Definition: GenEvent.cc:623
HepMC3::GenEvent::add_tree
void add_tree(const std::vector< GenParticlePtr > &particles)
Add whole tree in topological order.
Definition: GenEvent.cc:268
HepMC3::GenEvent::m_rootvertex
GenVertexPtr m_rootvertex
The root vertex is stored outside the normal vertices list to block user access to it.
Definition: GenEvent.h:359
HepMC3::GenEvent::length_unit
const Units::LengthUnit & length_unit() const
Get length unit.
Definition: GenEvent.h:143
HepMC3::GenEvent::attribute_as_string
string attribute_as_string(const string &name, const int &id=0) const
Get attribute of any type as string.
Definition: GenEvent.cc:800
HepMC3::GenEvent::remove_particles
void remove_particles(std::vector< GenParticlePtr > v)
Remove a set of particles.
Definition: GenEvent.cc:186
Errors.h
Implementation of error and warning macros.
GenRunInfo.h
Definition of class GenRunInfo.
HepMC3::GenEvent::m_particles
std::vector< GenParticlePtr > m_particles
List of particles.
Definition: GenEvent.h:343
HepMC3::GenEvent::attribute
shared_ptr< T > attribute(const string &name, const int &id=0) const
Get attribute of type T.
Definition: GenEvent.h:388
HepMC3::GenEvent::weight_names
const std::vector< std::string > & weight_names(const std::string &) const
Definition: GenEvent.h:111
HepMC3::GenVertex
Stores vertex-related information.
Definition: GenVertex.h:27
HepMC3::GenEvent::GenEvent
GenEvent(Units::MomentumUnit momentum_unit=Units::GEV, Units::LengthUnit length_unit=Units::MM)
Event constructor without a run.
Definition: GenEvent.cc:22
HepMC3::GenEvent::att_val_t
std::map< int, shared_ptr< Attribute > >::value_type att_val_t
Attribute map value type.
Definition: GenEvent.h:373
HepMC3::GenEvent::remove_vertex
void remove_vertex(GenVertexPtr v)
Remove vertex from the event.
Definition: GenEvent.cc:195
Units.h
Definition of class Units.
HepMC3::GenEvent::set_cross_section
void set_cross_section(GenCrossSectionPtr cs)
Set cross-section information.
Definition: GenEvent.h:167
HepMC3::GenEvent::weight
double weight(const std::string &name) const
Definition: GenEvent.h:96
HepMC3::GenEvent::~GenEvent
~GenEvent()
Destructor.
Definition: GenEvent.cc:76
HepMC3::GenEvent::m_length_unit
Units::LengthUnit m_length_unit
Length unit.
Definition: GenEvent.h:356
HepMC3::GenEvent::add_particle
void add_particle(GenParticlePtr p)
Add particle.
Definition: GenEvent.cc:50
HepMC3::GenEvent::remove_particle
void remove_particle(GenParticlePtr v)
Remove particle from the event.
Definition: GenEvent.cc:118
HepMC3::GenEvent::set_units
void set_units(Units::MomentumUnit new_momentum_unit, Units::LengthUnit new_length_unit)
Change event units Converts event from current units to new ones.
Definition: GenEvent.cc:396
HepMC3::GenParticle
Stores particle-related information.
Definition: GenParticle.h:31
HepMC3::GenEvent::read_data
void read_data(const GenEventData &data)
Fill GenEvent based on GenEventData.
Definition: GenEvent.cc:705
HepMC3::GenEvent::particles
const std::vector< ConstGenParticlePtr > & particles() const
Get list of particles (const)
Definition: GenEvent.cc:40
HepMC3::GenEvent::m_vertices
std::vector< GenVertexPtr > m_vertices
List of vertices.
Definition: GenEvent.h:345
HepMC3::GenEvent::att_key_t
std::map< string, std::map< int, shared_ptr< Attribute > > >::value_type att_key_t
Attribute map key type.
Definition: GenEvent.h:370
HepMC3::GenEvent::m_attributes
std::map< string, std::map< int, shared_ptr< Attribute > > > m_attributes
Map of event, particle and vertex attributes.
Definition: GenEvent.h:367
HepMC3::Units::LengthUnit
LengthUnit
Position units.
Definition: Units.h:32
GenCrossSection.h
Definition of attribute class GenCrossSection.
HepMC3::GenEvent::set_heavy_ion
void set_heavy_ion(GenHeavyIonPtr hi)
Set heavy ion generator additional information.
Definition: GenEvent.h:153
HepMC3::GenEvent::weight
double weight(const size_t &index=0) const
Definition: GenEvent.h:92
HepMC3::GenEvent::attribute_names
std::vector< string > attribute_names(const int &id=0) const
Get list of attribute names.
Definition: GenEvent.cc:635