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