HepMC3 event record library
HepMCCompatibility.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/// @file HepMCCompatibility.h
7/// @brief Implementation of compatibility layer (in-memory conversion functions) between HepMC2 and HepMC3
8///
9#ifndef HEPMCCOMPATIBILITY_H
10#define HEPMCCOMPATIBILITY_H
11
12#include "HepMC3/GenVertex.h"
13#include "HepMC3/GenParticle.h"
14#include "HepMC3/GenEvent.h"
15
16///Please note the HEPMC_HAS_CENTRALITY should be defined externaly
17#include "HepMC/GenVertex.h"
18#include "HepMC/GenParticle.h"
19#include "HepMC/GenEvent.h"
20
21/** Converts HepMC3::Genevent to HepMC::Genevent */
22HepMC::GenEvent* ConvertHepMCGenEvent_3to2(const HepMC3::GenEvent& evt)
23{
24 HepMC::GenEvent* n=new HepMC::GenEvent();
25 HepMC::Units::LengthUnit lu;
26 HepMC::Units::MomentumUnit mu;
27 if (evt.length_unit()==HepMC3::Units::CM) lu=HepMC::Units::CM;
28 if (evt.length_unit()==HepMC3::Units::MM) lu=HepMC::Units::MM;
29 if (evt.momentum_unit()==HepMC3::Units::MEV) mu=HepMC::Units::MEV;
30 if (evt.momentum_unit()==HepMC3::Units::GEV) mu=HepMC::Units::GEV;
31 n->use_units(mu,lu);
32
33
34 std::shared_ptr<HepMC3::IntAttribute> A_signal_process_vertex=evt.attribute<HepMC3::IntAttribute>("signal_process_vertex");
35 int signal_process_vertex_id=0;
36 if (A_signal_process_vertex) signal_process_vertex_id=A_signal_process_vertex->value();
37
38 std::map<HepMC3::ConstGenVertexPtr,HepMC::GenVertex*> vertexmap3to2;
39 for (auto v3: evt.vertices())
40 {
41 HepMC3::FourVector pos3=v3->position();
42 HepMC::FourVector pos2=HepMC::FourVector(pos3.x(),pos3.y(),pos3.pz(),pos3.t());
43 HepMC::GenVertex* v2= new HepMC::GenVertex(pos2);
44 std::vector<double> vweights;
45 if(v3->attribute_names().size())
46 {
47 std::shared_ptr<HepMC3::VectorDoubleAttribute> rsvec=v3->attribute<HepMC3::VectorDoubleAttribute>("weights");
48 if (rsvec) { vweights =rsvec->value(); }
49 else {
50 for (size_t ii=0; ii<100; ii++)
51 {
52 std::shared_ptr<HepMC3::DoubleAttribute> rs=v3->attribute<HepMC3::DoubleAttribute>("weight"+std::to_string((long long unsigned int)ii));
53 if (!rs) break;
54 vweights.push_back(rs->value());
55 }
56 }
57 }
58 if (vweights.empty())vweights.push_back(1.0);
59 v2->weights()=vweights;
60 v2->suggest_barcode(v3->id());
61 n->add_vertex(v2);
62 if (signal_process_vertex_id==v3->id()) n->set_signal_process_vertex(v2);
63 vertexmap3to2[v3]=v2;
64 }
65 std::map<HepMC3::ConstGenParticlePtr,HepMC::GenParticle*> particlemap3to2;
66 for (auto p3: evt.particles())
67 {
68 HepMC3::FourVector mom3=p3->momentum();
69 HepMC::FourVector mom2=HepMC::FourVector(mom3.px(),mom3.py(),mom3.pz(),mom3.e());
70 HepMC::GenParticle* p2= new HepMC::GenParticle(mom2,p3->pid(),p3->status());
71 if (p3->is_generated_mass_set())p2->set_generated_mass(p3->generated_mass());
72 /** Converts HepMC3::Genevent to HepMC::Genevent */
73 p2->suggest_barcode(10000+p3->id());
74 particlemap3to2[p3]=p2;
75
76 auto v3production=p3->production_vertex();
77 HepMC::GenVertex* v2production=nullptr;
78 if (v3production)
79 {
80 auto v2=vertexmap3to2.find(v3production);
81 if (v2!=vertexmap3to2.end()) v2production=vertexmap3to2[v3production];
82 }
83
84 auto v3end=p3->end_vertex();
85 HepMC::GenVertex* v2end=nullptr;
86 if (v3end)
87 {
88 auto v2=vertexmap3to2.find(v3end);
89 if (v2!=vertexmap3to2.end()) v2end=vertexmap3to2[v3end];
90 }
91
92 if (v2production) v2production->add_particle_out(p2);
93 if (v2end) v2end->add_particle_in(p2);
94
95 std::shared_ptr<HepMC3::DoubleAttribute> p3_theta=p3->attribute<HepMC3::DoubleAttribute>("theta");
96 std::shared_ptr<HepMC3::DoubleAttribute> p3_phi=p3->attribute<HepMC3::DoubleAttribute>("phi");
97 if (p3_phi && p3_theta)p2->set_polarization(HepMC::Polarization(p3_theta->value(),p3_phi->value()));
98 std::shared_ptr<HepMC3::VectorIntAttribute> flows = p3->attribute<HepMC3::VectorIntAttribute>("flows");
99 if (flows)
100 {
101 std::vector<int> flowvalues=flows->value();
102 for (size_t i=0; i<flowvalues.size(); i++) p2->set_flow(i+1,flowvalues.at(i));
103 } else {
104 std::shared_ptr<HepMC3::IntAttribute> p3_flow1=p3->attribute<HepMC3::IntAttribute>("flow1");
105 std::shared_ptr<HepMC3::IntAttribute> p3_flow2=p3->attribute<HepMC3::IntAttribute>("flow2");
106 std::shared_ptr<HepMC3::IntAttribute> p3_flow3=p3->attribute<HepMC3::IntAttribute>("flow3");
107 if (p3_flow1) p2->set_flow(1,p3_flow1->value());
108 if (p3_flow2) p2->set_flow(2,p3_flow2->value());
109 if (p3_flow3) p2->set_flow(3,p3_flow2->value());
110 }
111 }
112 std::vector<HepMC3::ConstGenParticlePtr> bms=evt.beams();
113 if (bms.size()>=2)
114 {
115 if (particlemap3to2.find(bms[0])!=particlemap3to2.end() && particlemap3to2.find(bms[1])!=particlemap3to2.end() )
116 {
117 n->set_beam_particles(particlemap3to2[bms[0]],particlemap3to2[bms[1]]);
118 }
119 }
120 std::shared_ptr<HepMC3::GenCrossSection> cs3 = evt.attribute<HepMC3::GenCrossSection>("GenCrossSection");
121 if(cs3) {
122 HepMC::GenCrossSection cs2;
123 cs2.set_cross_section(cs3->xsec(),cs3->xsec_err());
124 n->set_cross_section(cs2);
125 }
126 std::shared_ptr<HepMC3::GenPdfInfo> pdf3 = evt.attribute<HepMC3::GenPdfInfo>("GenPdfInfo");
127 if(pdf3)
128 {
129 HepMC::PdfInfo pdf2(
130 pdf3->parton_id[0],
131 pdf3->parton_id[1],
132 pdf3->x[0],
133 pdf3->x[1],
134 pdf3->scale,
135 pdf3->xf[0],
136 pdf3->xf[1],
137 pdf3->pdf_id[0],
138 pdf3->pdf_id[1]
139 );
140 n->set_pdf_info(pdf2);
141 }
142
143 std::shared_ptr<HepMC3::GenHeavyIon> hi3 = evt.attribute<HepMC3::GenHeavyIon>("GenHeavyIon");
144 if(hi3)
145 {
146 HepMC::HeavyIon hi2;
147 hi2.set_Ncoll_hard(hi3->Ncoll_hard);
148 hi2.set_Npart_proj(hi3->Npart_proj);
149 hi2.set_Npart_targ(hi3->Npart_targ);
150 hi2.set_Ncoll(hi3->Ncoll);
151 hi2.set_spectator_neutrons(hi3->spectator_neutrons);
152 hi2.set_spectator_protons(hi3->spectator_protons);
153 hi2.set_N_Nwounded_collisions(hi3->N_Nwounded_collisions);
154 hi2.set_Nwounded_N_collisions(hi3->Nwounded_N_collisions);
155 hi2.set_Nwounded_Nwounded_collisions(hi3->Nwounded_Nwounded_collisions);
156 hi2.set_impact_parameter(hi3->impact_parameter);
157 hi2.set_event_plane_angle(hi3->event_plane_angle);
158 hi2.set_eccentricity(hi3->eccentricity);
159 hi2.set_sigma_inel_NN(hi3->sigma_inel_NN);
160#ifdef HEPMC_HAS_CENTRALITY
161 hi2.set_centrality(hi3->centrality);
162#endif
163 n->set_heavy_ion(hi2);
164 }
165 std::vector<double> wv=evt.weights();
166 std::vector<std::string> wn=evt.weight_names();
167 for (size_t i=0; i<wv.size(); i++) n->weights()[wn.at(i)]=wv.at(i);
168
169 std::shared_ptr<HepMC3::DoubleAttribute> A_event_scale=evt.attribute<HepMC3::DoubleAttribute>("event_scale");
170 std::shared_ptr<HepMC3::DoubleAttribute> A_alphaQED=evt.attribute<HepMC3::DoubleAttribute>("alphaQED");
171 std::shared_ptr<HepMC3::DoubleAttribute> A_alphaQCD=evt.attribute<HepMC3::DoubleAttribute>("alphaQCD");
172 std::shared_ptr<HepMC3::IntAttribute> A_signal_process_id=evt.attribute<HepMC3::IntAttribute>("signal_process_id");
173 std::shared_ptr<HepMC3::IntAttribute> A_mpi=evt.attribute<HepMC3::IntAttribute>("mpi");
174
175
176 double event_scale=A_event_scale?(A_event_scale->value()):0.0;
177 double alphaQED=A_alphaQED?(A_alphaQED->value()):0.0;
178 double alphaQCD=A_alphaQCD?(A_alphaQCD->value()):0.0;
179 int signal_process_id=A_signal_process_id?(A_signal_process_id->value()):0;
180 int mpi=A_mpi?(A_mpi->value()):0;
181
182 std::vector<long> random_states;
183 for (int i=0; i<100; i++)
184 {
185 std::shared_ptr<HepMC3::IntAttribute> rs=evt.attribute<HepMC3::IntAttribute>("random_states"+std::to_string((long long unsigned int)i));
186 if (!rs) break;
187 random_states.push_back(rs->value());
188 }
189 n->set_mpi( mpi );
190 n->set_signal_process_id( signal_process_id );
191 n->set_event_number( evt.event_number() );
192 n->set_random_states( random_states );
193 n->set_event_scale( event_scale );
194 n->set_alphaQCD( alphaQCD );
195 n->set_alphaQED( alphaQED );
196 return n;
197}
198/** Converts HepMC::Genevent to HepMC3::Genevent */
199HepMC3::GenEvent* ConvertHepMCGenEvent_2to3(const HepMC::GenEvent& evt, std::shared_ptr<HepMC3::GenRunInfo> run )
200{
204 if (evt.length_unit()==HepMC::Units::CM) lu=HepMC3::Units::CM;
205 if (evt.length_unit()==HepMC::Units::MM) lu=HepMC3::Units::MM;
206 if (evt.momentum_unit()==HepMC::Units::MEV) mu=HepMC3::Units::MEV;
207 if (evt.momentum_unit()==HepMC::Units::GEV) mu=HepMC3::Units::GEV;
208 n->set_units(mu,lu);
209
210 n->set_run_info(run);
211 n->set_event_number(evt.event_number());
212 std::map<HepMC::GenVertex*,HepMC3::GenVertexPtr> vertexmap2to3;
213 for (auto v2=evt.vertices_begin(); v2!=evt.vertices_end(); v2++)
214 {
215 HepMC::FourVector pos2=(*v2)->position();
216 HepMC3::FourVector pos3=HepMC3::FourVector(pos2.x(),pos2.y(),pos2.pz(),pos2.t());
217 HepMC3::GenVertexPtr v3=std::make_shared<HepMC3::GenVertex>(pos3);
218 n->add_vertex(v3);
219 vertexmap2to3[*v2]=v3;
220 }
221
222 std::map<HepMC::GenParticle*,HepMC3::GenParticlePtr> particlemap2to3;
223 for (auto p2=evt.particles_begin(); p2!=evt.particles_end(); p2++)
224 {
225
226 HepMC::FourVector mom2=(*p2)->momentum();
227 HepMC3::FourVector mom3=HepMC3::FourVector(mom2.px(),mom2.py(),mom2.pz(),mom2.e());
228 HepMC3::GenParticlePtr p3= std::make_shared<HepMC3::GenParticle>(mom3,(*p2)->pdg_id(),(*p2)->status());
229 /// we set it always as there is no way to check if it is set
230 p3->set_generated_mass((*p2)->generated_mass());
231 particlemap2to3[*p2]=p3;
232
233 auto v2production=(*p2)->production_vertex();
234 HepMC3::GenVertexPtr v3production;
235 if (v2production)
236 {
237 auto v3=vertexmap2to3.find(v2production);
238 if (v3!=vertexmap2to3.end()) v3production=vertexmap2to3[v2production];
239 }
240
241 auto v2end=(*p2)->end_vertex();
242 HepMC3::GenVertexPtr v3end;
243 if (v2end)
244 {
245 auto v3=vertexmap2to3.find(v2end);
246 if (v3!=vertexmap2to3.end()) v3end=vertexmap2to3[v2end];
247 }
248
249 if (v3production) v3production->add_particle_out(p3);
250 if (v3end) v3end->add_particle_in(p3);
251
252 std::shared_ptr<HepMC3::DoubleAttribute> p3_theta=std::make_shared<HepMC3::DoubleAttribute>((*p2)->polarization().theta());
253
254 std::shared_ptr<HepMC3::DoubleAttribute> p3_phi=std::make_shared<HepMC3::DoubleAttribute>((*p2)->polarization().phi());
255 if (std::abs(p3_phi->value())>std::numeric_limits<double>::epsilon() || std::abs(p3_theta->value())>std::numeric_limits<double>::epsilon())
256 {
257 p3->add_attribute("theta",p3_theta);
258 p3->add_attribute("phi",p3_theta);
259 }
260
261 std::vector<int> flv;
262 auto fl=(*p2)->flow();
263 for (auto flel=fl.begin(); flel!=fl.end(); ++flel) flv.push_back((*flel).second);
264 std::shared_ptr<HepMC3::VectorIntAttribute> flows = std::make_shared<HepMC3::VectorIntAttribute>(flv);
265 if (flv.size())
266 p3->add_attribute("flows",flows);
267 }
268
269
270 std::shared_ptr<HepMC3::GenCrossSection> cs3 = std::make_shared<HepMC3::GenCrossSection>();
271 auto cs2=evt.cross_section();
272 if (cs2)cs3->set_cross_section(cs2->cross_section(),cs2->cross_section_error());
273 n->set_cross_section(cs3);
274
275
276 std::shared_ptr<HepMC3::GenPdfInfo> pdf3 = std::make_shared<HepMC3::GenPdfInfo>();
277 auto pdf2=evt.pdf_info();
278 if(pdf2)
279 {
280 pdf3->parton_id[0]=pdf2->id1();
281 pdf3->parton_id[1]=pdf2->id2();
282 pdf3->x[0]=pdf2->x1();
283 pdf3->x[1]=pdf2->x2();
284 pdf3->scale=pdf2->scalePDF();
285 pdf3->xf[0]=pdf2->pdf1();
286 pdf3->xf[1]=pdf2->pdf2();
287 pdf3->pdf_id[0]=pdf2->pdf_id1();
288 pdf3->pdf_id[1]=pdf2->pdf_id2();
289 n->set_pdf_info(pdf3);
290 }
291
292 std::shared_ptr<HepMC3::GenHeavyIon> hi3 = std::make_shared<HepMC3::GenHeavyIon>();
293 auto hi2=evt.heavy_ion();
294 if(hi2)
295 {
296 hi3->Ncoll_hard=hi2->Ncoll_hard();
297 hi3->Npart_proj=hi2->Npart_proj();
298 hi3->Npart_targ=hi2->Npart_targ();
299 hi3->Ncoll=hi2->Ncoll();
300 hi3->spectator_neutrons=hi2->spectator_neutrons();
301 hi3->spectator_protons=hi2->spectator_protons();
302 hi3->N_Nwounded_collisions=hi2->N_Nwounded_collisions();
303 hi3->Nwounded_N_collisions=hi2->Nwounded_N_collisions();
304 hi3->Nwounded_Nwounded_collisions=hi2->Nwounded_Nwounded_collisions();
305 hi3->impact_parameter=hi2->impact_parameter();
306 hi3->event_plane_angle=hi2->event_plane_angle();
307 hi3->eccentricity=hi2->eccentricity();
308 hi3->sigma_inel_NN=hi2->sigma_inel_NN();
309#ifdef HEPMC_HAS_CENTRALITY
310 hi3->centrality=hi2->centrality();
311#endif
312 n->set_heavy_ion(hi3);
313 }
314 /** Yes, the desing is not always perfect */
315 std::stringstream ss;
316 evt.weights().print(ss);
317 std::string allweights=ss.str();
318 std::vector<std::string> wnames;
319 std::string token;
320 std::size_t pos=0;
321 for (;;)
322 {
323 std::size_t name_begin=allweights.find_first_not_of(" (\n",pos);
324 if (name_begin==std::string::npos) break;
325 std::size_t name_end=allweights.find_first_of(",",name_begin);
326 wnames.push_back(allweights.substr(name_begin,name_end-name_begin));
327 pos=allweights.find_first_of(")",name_end)+1;
328 }
329 n->run_info()->set_weight_names(wnames);
330 n->weights()=std::vector<double>(wnames.size(),1.0);
331 for (auto wn: wnames)
332 {
333 n->weight(wn)=evt.weights()[wn];
334 }
335
336 std::shared_ptr<HepMC3::DoubleAttribute> A_event_scale=std::make_shared<HepMC3::DoubleAttribute>(evt.event_scale());
337 n->add_attribute("event_scale",A_event_scale);
338 std::shared_ptr<HepMC3::DoubleAttribute> A_alphaQED=std::make_shared<HepMC3::DoubleAttribute>(evt.alphaQED());
339 n->add_attribute("alphaQED",A_alphaQED);
340 std::shared_ptr<HepMC3::DoubleAttribute> A_alphaQCD=std::make_shared<HepMC3::DoubleAttribute>(evt.alphaQCD());
341 n->add_attribute("alphaQCD",A_alphaQCD);
342 std::shared_ptr<HepMC3::IntAttribute> A_signal_process_id=std::make_shared<HepMC3::IntAttribute>(evt.signal_process_id());
343 n->add_attribute("signal_process_id",A_signal_process_id);
344 std::shared_ptr<HepMC3::IntAttribute> A_mpi=std::make_shared<HepMC3::IntAttribute>(evt.mpi());
345 n->add_attribute("mpi",A_mpi);
346 std::shared_ptr<HepMC3::VectorLongIntAttribute> A_random_states= std::make_shared<HepMC3::VectorLongIntAttribute>( evt.random_states());
347 n->add_attribute("random_states",A_random_states);
348
349 auto spv2=evt.signal_process_vertex();
350 if (vertexmap2to3.find(spv2)!=vertexmap2to3.end())
351 {
352 int signal_process_vertex_id3=vertexmap2to3[spv2]->id();
353 std::shared_ptr<HepMC3::IntAttribute> A_signal_process_vertex=std::make_shared<HepMC3::IntAttribute>(signal_process_vertex_id3);
354 n->add_attribute("signal_process_vertex",A_signal_process_vertex);
355 }
356 return n;
357}
358#endif
Definition of class GenEvent.
Definition of class GenParticle.
Definition of class GenVertex.
HepMC::GenEvent * ConvertHepMCGenEvent_3to2(const HepMC3::GenEvent &evt)
Please note the HEPMC_HAS_CENTRALITY should be defined externaly.
HepMC3::GenEvent * ConvertHepMCGenEvent_2to3(const HepMC::GenEvent &evt, std::shared_ptr< HepMC3::GenRunInfo > run)
Attribute that holds a real number as a double.
Definition: Attribute.h:241
Generic 4-vector.
Definition: FourVector.h:36
double e() const
Energy component of momentum.
Definition: FourVector.h:131
double pz() const
z-component of momentum
Definition: FourVector.h:124
double t() const
Time component of position/displacement.
Definition: FourVector.h:102
double px() const
x-component of momentum
Definition: FourVector.h:110
double py() const
y-component of momentum
Definition: FourVector.h:117
double x() const
x-component of position/displacement
Definition: FourVector.h:81
double y() const
y-component of position/displacement
Definition: FourVector.h:88
Stores additional information about cross-section.
void set_cross_section(const double &xs, const double &xs_err, const long &n_acc=-1, const long &n_att=-1)
Set all fields.
Stores event-related information.
Definition: GenEvent.h:41
std::shared_ptr< T > attribute(const std::string &name, const int &id=0) const
Get attribute of type T.
Definition: GenEvent.h:409
void add_vertex(GenVertexPtr v)
Add vertex.
Definition: GenEvent.cc:96
int event_number() const
Get event number.
Definition: GenEvent.h:148
std::shared_ptr< GenRunInfo > run_info() const
Get a pointer to the the GenRunInfo object.
Definition: GenEvent.h:137
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:391
void set_event_number(const int &num)
Set event number.
Definition: GenEvent.h:150
std::vector< ConstGenParticlePtr > beams() const
Vector of beam particles.
Definition: GenEvent.cc:416
void set_run_info(std::shared_ptr< GenRunInfo > run)
Set the GenRunInfo object by smart pointer.
Definition: GenEvent.h:141
const std::vector< ConstGenVertexPtr > & vertices() const
Get list of vertices (const)
Definition: GenEvent.cc:43
const Units::MomentumUnit & momentum_unit() const
Get momentum unit.
Definition: GenEvent.h:153
const Units::LengthUnit & length_unit() const
Get length unit.
Definition: GenEvent.h:155
double weight(const unsigned long &index=0) const
Definition: GenEvent.h:103
void set_cross_section(GenCrossSectionPtr cs)
Set cross-section information.
Definition: GenEvent.h:179
void set_heavy_ion(GenHeavyIonPtr hi)
Set heavy ion generator additional information.
Definition: GenEvent.h:165
void add_attribute(const std::string &name, const std::shared_ptr< Attribute > &att, const int &id=0)
Definition: GenEvent.cc:806
const std::vector< std::string > & weight_names() const
Definition: GenEvent.h:123
const std::vector< double > & weights() const
Get event weight values as a vector.
Definition: GenEvent.h:98
void set_pdf_info(GenPdfInfoPtr pi)
Set PDF information.
Definition: GenEvent.h:172
const std::vector< ConstGenParticlePtr > & particles() const
Get list of particles (const)
Definition: GenEvent.cc:39
Stores additional information about Heavy Ion generator.
Definition: GenHeavyIon.h:28
Stores additional information about PDFs.
Definition: GenPdfInfo.h:32
Attribute that holds an Integer implemented as an int.
Definition: Attribute.h:157
int value() const
get the value associated to this Attribute.
Definition: Attribute.h:179
LengthUnit
Position units.
Definition: Units.h:32
MomentumUnit
Momentum units.
Definition: Units.h:29
Attribute that holds a vector of FPs of type double.
Definition: Attribute.h:1091
std::vector< double > value() const
get the value associated to this Attribute.
Definition: Attribute.h:1117
Attribute that holds a vector of integers of type int.
Definition: Attribute.h:1001
std::vector< int > value() const
get the value associated to this Attribute.
Definition: Attribute.h:1027