22 #include "HepPDT/ParticleDataTable.hh"
28 static const char *
const s_protocolSep =
":";
34 const std::string&
name,
35 const IInterface*
parent ) :
37 m_ioBackend( nullptr )
43 declareProperty(
"Output",
45 "Name of the back-end we'll use to write out the HepMC::GenEvent."
46 "\nEx: ascii:hepmc.genevent.txt" );
49 declareProperty(
"McEvents",
51 "Input location of the McEventCollection to write out" );
77 (*m_ioBackend) <<
"# epsilon<float> = "
78 << std::numeric_limits<float>::epsilon()
80 <<
"# dbl precision = "
81 << std::numeric_limits<double>::digits10 + 1
83 <<
"# flt precision = "
84 << std::numeric_limits<float>::digits10 + 1
87 return StatusCode::SUCCESS;
93 return StatusCode::SUCCESS;
102 return StatusCode::FAILURE;
105 if ( mcEvts->
empty() ) {
107 return StatusCode::FAILURE;
110 const HepMC::GenEvent *
evt = mcEvts->
front();
112 ATH_MSG_ERROR(
"Retrieved NULL pointer to HepMC::GenEvent !!");
113 return StatusCode::FAILURE;
125 std::ostringstream
out;
129 out.precision( std::numeric_limits<double>::digits10 + 1 );
131 out.setf(std::ios::dec,std::ios::basefield);
132 out.setf(std::ios::scientific,std::ios::floatfield);
136 long evt_vertices_size=
evt->vertices().size();
137 std::shared_ptr<HepMC3::DoubleAttribute> A_alphaQCD=
evt->attribute<HepMC3::DoubleAttribute>(
"alphaQCD");
138 double evt_alphaQCD=(A_alphaQCD?(A_alphaQCD->value()):0.0);
139 std::shared_ptr<HepMC3::DoubleAttribute> A_alphaQED=
evt->attribute<HepMC3::DoubleAttribute>(
"alphaQED");
140 double evt_alphaQED=(A_alphaQED?(A_alphaQED->value()):0.0);
141 std::shared_ptr<HepMC3::DoubleAttribute> A_event_scale=
evt->attribute<HepMC3::DoubleAttribute>(
"event_scale");
142 double evt_event_scale=(A_event_scale?(A_event_scale->value()):0.0);
143 std::shared_ptr<HepMC3::VectorLongIntAttribute> A_random_states=
evt->attribute<HepMC3::VectorLongIntAttribute>(
"random_states");
144 std::vector<long int> random_states=(A_random_states?(A_random_states->value()):std::vector<long int>());
145 long random_states_size=random_states.size();
147 long evt_vertices_size=
evt->vertices_size();
148 double evt_alphaQCD=
evt->alphaQCD();
149 double evt_alphaQED=
evt->alphaQED();
150 double evt_event_scale=
evt->event_scale();
151 std::vector<long int> random_states=
evt->random_states();
152 long random_states_size=random_states.size();
154 out <<
"# -- GenEvent -->\n";
155 out <<
"#" <<
evt->event_number()
156 <<
" " << evt_event_scale
157 <<
" " << evt_alphaQCD
158 <<
" " << evt_alphaQED
161 <<
" " << evt_vertices_size
162 <<
" " << random_states_size
165 std::copy( random_states.begin(), random_states.end(), std::ostream_iterator<long int>(
out,
" ") );
166 out <<
evt->weights().size() <<
"\n";
169 std::copy(
evt->weights().begin(),
evt->weights().end(), std::ostream_iterator<double>(
out,
" ") );
172 out <<
"#-- particles --\n";
173 for (
const auto&
p: *
evt) {
177 const HepMC::FourVector
mom =
p->momentum();
178 std::ostringstream buf;
179 buf.precision( std::numeric_limits<float>::digits10 + 1 );
180 buf.setf(std::ios::dec,std::ios::basefield);
181 buf.setf(std::ios::scientific,std::ios::floatfield);
183 const float px =
static_cast<float>(
mom.px());
184 const float py =
static_cast<float>(
mom.py());
185 const float pz =
static_cast<float>(
mom.pz());
186 const float m =
static_cast<float>(
mom.m());
188 buf <<
px <<
" " <<
py <<
" " <<
pz <<
" " <<
e <<
" " <<
m <<
"\n";
192 out <<
"# "<<
p->status()
193 <<
" " << pol.theta()
201 out <<
"#-- vertices -- \n";
203 for (
const auto&
v:
evt->vertices()) {
206 const HepMC::FourVector
pos =
v->position();
207 std::ostringstream buf;
208 buf.precision( std::numeric_limits<float>::digits10 + 1 );
209 buf.setf(std::ios::dec,std::ios::basefield);
210 buf.setf(std::ios::scientific,std::ios::floatfield);
212 buf <<
pos.x() <<
" " <<
pos.y() <<
" " <<
pos.z() <<
" " <<
pos.t() <<
"\n";
216 std::string svertexeights(
"1.0");
217 auto vertexeights=
v->attribute<HepMC3::VectorDoubleAttribute>(
"weights");
218 if (vertexeights) vertexeights->to_string(svertexeights);
219 out << svertexeights;
223 for ( HepMC::GenEvent::vertex_const_iterator
224 i =
evt->vertices_begin(),
225 iEnd =
evt->vertices_end();
228 const HepMC::GenVertex *
v = *
i;
230 out <<
"# " <<
v->barcode() <<
" " <<
v->id() <<
"\n";
232 const HepMC::FourVector
pos =
v->position();
233 std::ostringstream buf;
234 buf.precision( std::numeric_limits<float>::digits10 + 1 );
235 buf.setf(std::ios::dec,std::ios::basefield);
236 buf.setf(std::ios::scientific,std::ios::floatfield);
238 buf <<
pos.x() <<
" " <<
pos.y() <<
" " <<
pos.z() <<
" " <<
pos.t() <<
"\n";
242 std::copy(
v->weights().begin(),
v->weights().end(), std::ostream_iterator<double>(
out,
" ") );
247 out <<
"#<-- GenEvent --\n";
250 return StatusCode::SUCCESS;
260 std::string protocol =
"ascii";
261 std::string
fileName =
"hepmc.genevent.txt";
272 std::string::size_type protocolPos =
url.find(s_protocolSep);
274 if ( std::string::npos != protocolPos ) {
275 protocol =
url.substr( 0, protocolPos );
276 fileName =
url.substr( protocolPos+1, std::string::npos );
283 std::transform( protocol.begin(), protocol.end(), protocol.begin(), [](
unsigned char c){ return std::tolower(c); } );
284 if (
"ascii" == protocol ) {
288 ATH_MSG_WARNING(
"UNKNOWN protocol [" << protocol <<
"] !!" <<
endmsg <<
"Will use [ascii] instead...");