HepMC3 event record library
HEPEVT_Wrapper.h
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 #ifndef HEPMC3_HEPEVT_WRAPPER_H
7 #define HEPMC3_HEPEVT_WRAPPER_H
8 /**
9  * @file HEPEVT_Wrapper.h
10  * @brief Definition of \b class HEPEVT_Wrapper
11  *
12  * @class HepMC3::HEPEVT_Wrapper
13  * @brief An interface to HEPEVT common block implemented in a traditional way.
14  * When possible this implementation should be avoided and the templated
15  * version should be used instead.
16  *
17  * @note This header file does not include HEPEVT definition, only declaration.
18  * Including this wrapper requires that HEPEVT is defined somewhere
19  * in the project (most likely as FORTRAN common block).
20  *
21  * @note Make sure that HEPEVT definition in project matches this definition
22  * (NMXHEP, double precision, etc.) Change this definition if necessary.
23  *
24  */
25 
26 #ifndef HEPMC3_HEPEVT_NMXHEP
27 /** Default number of particles in the HEPEVT structure */
28 #define HEPMC3_HEPEVT_NMXHEP 10000
29 #endif
30 
31 #ifndef HEPMC3_HEPEVT_PRECISION
32 /** Default precision of the 4-momentum, time-space position and mass */
33 #define HEPMC3_HEPEVT_PRECISION double
34 #endif
35 
36 /* This definition of HEPEVT corresponds to FORTRAN definition:
37 
38  PARAMETER (NMXHEP=10000)
39  COMMON /HEPEVT/ NEVHEP,NHEP,ISTHEP(NMXHEP),IDHEP(NMXHEP),
40  & JMOHEP(2,NMXHEP),JDAHEP(2,NMXHEP),PHEP(5,NMXHEP),VHEP(4,NMXHEP)
41  INTEGER NEVHEP,NHEP,ISTHEP,IDHEP,JMOHEP,JDAHEP
42  DOUBLE PRECISION PHEP,VHEP
43 */
44 
45 static const int NMXHEP = HEPMC3_HEPEVT_NMXHEP; //!< Number of particles in the HEPEVT structure
46 typedef HEPMC3_HEPEVT_PRECISION momentum_t; //!< Precision of the 4-momentum, time-space position and mass
47 
48 /** @struct HEPEVT
49  * @brief C structure representing Fortran common block HEPEVT
50  * T. Sjöstrand et al., "A proposed standard event record",
51  * in `Z physics at LEP 1', eds. G. Altarelli, R. Kleiss and C. Verzegnassi,
52  * Geneva, Switzerland, September 4-5, 1989, CERN 89-08 (Geneva, 1989), Vol. 3, p. 327
53  * Disk representation is given by Fortran WRITE/READ format.
54  */
55 struct HEPEVT
56 {
57  int nevhep; //!< Event number
58  int nhep; //!< Number of entries in the event
59  int isthep[NMXHEP]; //!< Status code
60  int idhep [NMXHEP]; //!< PDG ID
61  int jmohep[NMXHEP][2]; //!< Pointer to position of 1st and 2nd (or last!) mother
62  int jdahep[NMXHEP][2]; //!< Pointer to position of 1nd and 2nd (or last!) daughter
63  momentum_t phep [NMXHEP][5]; //!< Momentum: px, py, pz, e, m
64  momentum_t vhep [NMXHEP][4]; //!< Time-space position: x, y, z, t
65 }; //!< Fortran common block HEPEVT
66 
67 #include <iostream>
68 #include <cstdio>
69 #include <set>
70 #include <map>
71 #include <cstring> // memset
72 #include <cassert>
73 #include <algorithm> //min max for VS2017
74 #ifndef HEPEVT_WRAPPER_HEADER_ONLY
75 #include "HepMC3/GenEvent.h"
76 #include "HepMC3/GenParticle.h"
77 #include "HepMC3/GenVertex.h"
78 #include "HepMC3/HEPEVT_Helpers.h"
79 #endif
80 /** @brief Pointer to HEPEVT common block */
81 HEPMC3_EXPORT_API struct HEPEVT* hepevtptr;
82 
83 namespace HepMC3
84 {
86 {
87 
88 //
89 // Functions
90 //
91 public:
92 
93  /** @brief Print information from HEPEVT common block */
94  static void print_hepevt( std::ostream& ostr = std::cout );
95  /** @brief Print particle information */
96  static void print_hepevt_particle( int index, std::ostream& ostr = std::cout );
97  /** @brief Set all entries in HEPEVT to zero */
98  static void zero_everything();
99 #ifndef HEPEVT_WRAPPER_HEADER_ONLY
100  /** @brief Convert GenEvent to HEPEVT*/
101  static bool GenEvent_to_HEPEVT( const GenEvent* evt ) { return GenEvent_to_HEPEVT_static<HEPEVT_Wrapper>(evt);};
102  /** @brief Convert HEPEVT to GenEvent*/
103  static bool HEPEVT_to_GenEvent( GenEvent* evt ) { return HEPEVT_to_GenEvent_static<HEPEVT_Wrapper>(evt);};
104 #endif
105  /** @brief Tries to fix list of daughters */
106  static bool fix_daughters();
107 //
108 // Accessors
109 //
110 public:
111  static void set_max_number_entries( unsigned int size ) { if (size != NMXHEP) printf("This implementation does not support change of the block size.\n"); assert(size == NMXHEP); }//!< Set block size
112  static void set_hepevt_address(char *c) { hepevtptr = (struct HEPEVT*)c; } //!< Set Fortran block address
113  static int max_number_entries() { return NMXHEP; } //!< Block size
114  static int event_number() { return hepevtptr->nevhep; } //!< Get event number
115  static int number_entries() { return hepevtptr->nhep; } //!< Get number of entries
116  static int status(const int& index ) { return hepevtptr->isthep[index-1]; } //!< Get status code
117  static int id(const int& index ) { return hepevtptr->idhep[index-1]; } //!< Get PDG particle id
118  static int first_parent(const int& index ) { return hepevtptr->jmohep[index-1][0]; } //!< Get index of 1st mother
119  static int last_parent(const int& index ) { return hepevtptr->jmohep[index-1][1]; } //!< Get index of last mother
120  static int first_child(const int& index ) { return hepevtptr->jdahep[index-1][0]; } //!< Get index of 1st daughter
121  static int last_child(const int& index ) { return hepevtptr->jdahep[index-1][1]; } //!< Get index of last daughter
122  static double px(const int& index ) { return hepevtptr->phep[index-1][0]; } //!< Get X momentum
123  static double py(const int& index ) { return hepevtptr->phep[index-1][1]; } //!< Get Y momentum
124  static double pz(const int& index ) { return hepevtptr->phep[index-1][2]; } //!< Get Z momentum
125  static double e(const int& index ) { return hepevtptr->phep[index-1][3]; } //!< Get Energy
126  static double m(const int& index ) { return hepevtptr->phep[index-1][4]; } //!< Get generated mass
127  static double x(const int& index ) { return hepevtptr->vhep[index-1][0]; } //!< Get X Production vertex
128  static double y(const int& index ) { return hepevtptr->vhep[index-1][1]; } //!< Get Y Production vertex
129  static double z(const int& index ) { return hepevtptr->vhep[index-1][2]; } //!< Get Z Production vertex
130  static double t(const int& index ) { return hepevtptr->vhep[index-1][3]; } //!< Get production time
131  static int number_parents(const int& index ); //!< Get number of parents
132  static int number_children(const int& index ); //!< Get number of children from the range of daughters
133  static int number_children_exact(const int& index ); //!< Get number of children by counting
134  static void set_event_number( const int& evtno ) { hepevtptr->nevhep = evtno; } //!< Set event number
135  static void set_number_entries( const int& noentries ) { hepevtptr->nhep = noentries; } //!< Set number of entries
136  static void set_status( const int& index, const int& status ) { hepevtptr->isthep[index-1] = status; } //!< Set status code
137  static void set_id(const int& index, const int& id ) { hepevtptr->idhep[index-1] = id; } //!< Set PDG particle id
138  static void set_parents( const int& index, const int& firstparent, const int& lastparent ); //!< Set parents
139  static void set_children( const int& index, const int& firstchild, const int& lastchild ); //!< Set children
140  static void set_momentum( const int& index, const double& px, const double& py, const double& pz, const double& e ); //!< Set 4-momentum
141  static void set_mass( const int& index, double mass ); //!< Set mass
142  static void set_position( const int& index, const double& x, const double& y, const double& z, const double& t ); //!< Set position in time-space
143 };
144 
145 //
146 // inline definitions
147 //
148 inline void HEPEVT_Wrapper::print_hepevt( std::ostream& ostr )
149 {
150  ostr << " Event No.: " << hepevtptr->nevhep << std::endl;
151  ostr<< " Nr Type Parent(s) Daughter(s) Px Py Pz E Inv. M." << std::endl;
152  for ( int i = 1; i <= hepevtptr->nhep; ++i )
153  {
155  }
156 }
157 
158 inline void HEPEVT_Wrapper::print_hepevt_particle( int index, std::ostream& ostr )
159 {
160  char buf[255];//Note: the format is fixed, so no reason for complicated treatment
161 
162  sprintf(buf,"%5i %6i",index,hepevtptr->idhep[index-1]);
163  ostr << buf;
164  sprintf(buf,"%4i - %4i ",hepevtptr->jmohep[index-1][0],hepevtptr->jmohep[index-1][1]);
165  ostr << buf;
166  sprintf(buf,"%4i - %4i ",hepevtptr->jdahep[index-1][0],hepevtptr->jdahep[index-1][1]);
167  ostr << buf;
168  sprintf(buf,"%8.2f %8.2f %8.2f %8.2f %8.2f",hepevtptr->phep[index-1][0],hepevtptr->phep[index-1][1],hepevtptr->phep[index-1][2],hepevtptr->phep[index-1][3],hepevtptr->phep[index-1][4]);
169  ostr << buf << std::endl;
170 }
171 
173 {
174  memset(hepevtptr,0,sizeof(struct HEPEVT));
175 }
176 
177 inline int HEPEVT_Wrapper::number_parents( const int& index )
178 {
179  return (hepevtptr->jmohep[index-1][0]) ? (hepevtptr->jmohep[index-1][1]) ? hepevtptr->jmohep[index-1][1]-hepevtptr->jmohep[index-1][0] : 1 : 0;
180 }
181 
182 inline int HEPEVT_Wrapper::number_children( const int& index )
183 {
184  return (hepevtptr->jdahep[index-1][0]) ? (hepevtptr->jdahep[index-1][1]) ? hepevtptr->jdahep[index-1][1]-hepevtptr->jdahep[index-1][0] : 1 : 0;
185 }
186 
187 inline int HEPEVT_Wrapper::number_children_exact( const int& index )
188 {
189  int nc=0;
190  for ( int i = 1; i <= hepevtptr->nhep; ++i )
191  if (((hepevtptr->jmohep[i-1][0] <= index && hepevtptr->jmohep[i-1][1]>=index))||(hepevtptr->jmohep[i-1][0]==index)||(hepevtptr->jmohep[i-1][1]==index)) nc++;
192  return nc;
193 }
194 
195 inline void HEPEVT_Wrapper::set_parents( const int& index, const int& firstparent, const int&lastparent )
196 {
197  hepevtptr->jmohep[index-1][0] = firstparent;
198  hepevtptr->jmohep[index-1][1] = lastparent;
199 }
200 
201 inline void HEPEVT_Wrapper::set_children( const int& index, const int& firstchild, const int& lastchild )
202 {
203  hepevtptr->jdahep[index-1][0] = firstchild;
204  hepevtptr->jdahep[index-1][1] = lastchild;
205 }
206 
207 inline void HEPEVT_Wrapper::set_momentum( const int& index, const double& px,const double& py, const double& pz, const double& e )
208 {
209  hepevtptr->phep[index-1][0] = px;
210  hepevtptr->phep[index-1][1] = py;
211  hepevtptr->phep[index-1][2] = pz;
212  hepevtptr->phep[index-1][3] = e;
213 }
214 
215 inline void HEPEVT_Wrapper::set_mass( const int& index, double mass )
216 {
217  hepevtptr->phep[index-1][4] = mass;
218 }
219 
220 inline void HEPEVT_Wrapper::set_position( const int& index, const double& x, const double& y, const double& z, const double& t )
221 {
222  hepevtptr->vhep[index-1][0] = x;
223  hepevtptr->vhep[index-1][1] = y;
224  hepevtptr->vhep[index-1][2] = z;
225  hepevtptr->vhep[index-1][3] = t;
226 }
227 
229 {
230  /*AV The function should be called for a record that has correct particle ordering and mother ids.
231  As a result it produces a record with ranges where the daughters can be found.
232  Not every particle in the range will be a daughter. It is true only for proper events.
233  The return tells if the record was fixed succesfully.
234  */
235  for ( int i = 1; i <= HEPEVT_Wrapper::number_entries(); i++ )
236  for ( int k = 1; k <= HEPEVT_Wrapper::number_entries(); k++ ) if (i!=k)
239  bool is_fixed = true;
240  for ( int i = 1; i <= HEPEVT_Wrapper::number_entries(); i++ )
242  return is_fixed;
243 }
244 
245 } // namespace HepMC3
246 #endif
Definition of class GenEvent.
Definition of class GenParticle.
Definition of class GenVertex.
Helper functions used to manipulate with HEPEVT block.
#define HEPMC3_HEPEVT_NMXHEP
HEPMC3_HEPEVT_PRECISION momentum_t
Precision of the 4-momentum, time-space position and mass.
#define HEPMC3_HEPEVT_PRECISION
static const int NMXHEP
Number of particles in the HEPEVT structure.
HEPMC3_EXPORT_API struct HEPEVT * hepevtptr
Pointer to HEPEVT common block.
Stores event-related information.
Definition: GenEvent.h:41
An interface to HEPEVT common block implemented in a traditional way. When possible this implementati...
static double pz(const int &index)
Get Z momentum.
static int last_child(const int &index)
Get index of last daughter.
static void set_mass(const int &index, double mass)
Set mass.
static double py(const int &index)
Get Y momentum.
static void set_momentum(const int &index, const double &px, const double &py, const double &pz, const double &e)
Set 4-momentum.
static void set_number_entries(const int &noentries)
Set number of entries.
static int number_children_exact(const int &index)
Get number of children by counting.
static double m(const int &index)
Get generated mass.
static int number_parents(const int &index)
Get number of parents.
static bool HEPEVT_to_GenEvent(GenEvent *evt)
Convert HEPEVT to GenEvent.
static void print_hepevt_particle(int index, std::ostream &ostr=std::cout)
Print particle information.
static void print_hepevt(std::ostream &ostr=std::cout)
Print information from HEPEVT common block.
static void zero_everything()
Set all entries in HEPEVT to zero.
static void set_max_number_entries(unsigned int size)
Set block size.
static void set_id(const int &index, const int &id)
Set PDG particle id.
static int max_number_entries()
Block size.
static void set_status(const int &index, const int &status)
Set status code.
static bool fix_daughters()
Tries to fix list of daughters.
static void set_hepevt_address(char *c)
Set Fortran block address.
static double y(const int &index)
Get Y Production vertex.
static double t(const int &index)
Get production time.
static int last_parent(const int &index)
Get index of last mother.
static double e(const int &index)
Get Energy.
static double z(const int &index)
Get Z Production vertex.
static double x(const int &index)
Get X Production vertex.
static bool GenEvent_to_HEPEVT(const GenEvent *evt)
Convert GenEvent to HEPEVT.
static int first_child(const int &index)
Get index of 1st daughter.
static int event_number()
Get event number.
static int status(const int &index)
Get status code.
static void set_parents(const int &index, const int &firstparent, const int &lastparent)
Set parents.
static void set_position(const int &index, const double &x, const double &y, const double &z, const double &t)
Set position in time-space.
static double px(const int &index)
Get X momentum.
static int first_parent(const int &index)
Get index of 1st mother.
static void set_children(const int &index, const int &firstchild, const int &lastchild)
Set children.
static int number_entries()
Get number of entries.
static int id(const int &index)
Get PDG particle id.
static void set_event_number(const int &evtno)
Set event number.
static int number_children(const int &index)
Get number of children from the range of daughters.
HepMC3 main namespace.
Fortran common block HEPEVT.
int jmohep[NMXHEP][2]
Pointer to position of 1st and 2nd (or last!) mother.
int isthep[NMXHEP]
Status code.
int idhep[NMXHEP]
PDG ID.
momentum_t phep[NMXHEP][5]
Momentum: px, py, pz, e, m.
momentum_t vhep[NMXHEP][4]
Time-space position: x, y, z, t.
int jdahep[NMXHEP][2]
Pointer to position of 1nd and 2nd (or last!) daughter.
int nevhep
Event number.
int nhep
Number of entries in the event.