9 #ifndef HEPMCCOMPATIBILITY_H
10 #define HEPMCCOMPATIBILITY_H
17 #include "HepMC/GenVertex.h"
18 #include "HepMC/GenParticle.h"
19 #include "HepMC/GenEvent.h"
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;
35 int signal_process_vertex_id=0;
36 if (A_signal_process_vertex) signal_process_vertex_id=A_signal_process_vertex->
value();
38 std::map<HepMC3::ConstGenVertexPtr,HepMC::GenVertex*> vertexmap3to2;
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())
48 if (rsvec) { vweights =rsvec->value(); }
50 for (
size_t ii=0; ii<100; ii++)
52 std::shared_ptr<HepMC3::DoubleAttribute> rs=v3->attribute<
HepMC3::DoubleAttribute>(
"weight"+std::to_string((
long long unsigned int)ii));
54 vweights.push_back(rs->value());
58 if (vweights.empty())vweights.push_back(1.0);
59 v2->weights()=vweights;
60 v2->suggest_barcode(v3->id());
62 if (signal_process_vertex_id==v3->id()) n->set_signal_process_vertex(v2);
65 std::map<HepMC3::ConstGenParticlePtr,HepMC::GenParticle*> particlemap3to2;
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());
73 p2->suggest_barcode(10000+p3->id());
74 particlemap3to2[p3]=p2;
76 auto v3production=p3->production_vertex();
77 HepMC::GenVertex* v2production=
nullptr;
80 auto v2=vertexmap3to2.find(v3production);
81 if (v2!=vertexmap3to2.end()) v2production=vertexmap3to2[v3production];
84 auto v3end=p3->end_vertex();
85 HepMC::GenVertex* v2end=
nullptr;
88 auto v2=vertexmap3to2.find(v3end);
89 if (v2!=vertexmap3to2.end()) v2end=vertexmap3to2[v3end];
92 if (v2production) v2production->add_particle_out(p2);
93 if (v2end) v2end->add_particle_in(p2);
97 if (p3_phi && p3_theta)p2->set_polarization(HepMC::Polarization(p3_theta->value(),p3_phi->value()));
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));
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());
112 std::vector<HepMC3::ConstGenParticlePtr> bms=evt.
beams();
115 if (particlemap3to2.find(bms[0])!=particlemap3to2.end() && particlemap3to2.find(bms[1])!=particlemap3to2.end() )
117 n->set_beam_particles(particlemap3to2[bms[0]],particlemap3to2[bms[1]]);
122 HepMC::GenCrossSection cs2;
124 n->set_cross_section(cs2);
140 n->set_pdf_info(pdf2);
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);
163 n->set_heavy_ion(hi2);
165 std::vector<double> wv=evt.
weights();
167 for (
size_t i=0; i<wv.size(); i++) n->weights()[wn.at(i)]=wv.at(i);
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;
182 std::vector<long> random_states;
183 for (
int i=0; i<100; i++)
187 random_states.push_back(rs->value());
190 n->set_signal_process_id( signal_process_id );
192 n->set_random_states( random_states );
193 n->set_event_scale( event_scale );
194 n->set_alphaQCD( alphaQCD );
195 n->set_alphaQED( alphaQED );
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;
212 std::map<HepMC::GenVertex*,HepMC3::GenVertexPtr> vertexmap2to3;
213 for (
auto v2=evt.vertices_begin(); v2!=evt.vertices_end(); v2++)
215 HepMC::FourVector pos2=(*v2)->position();
217 HepMC3::GenVertexPtr v3=std::make_shared<HepMC3::GenVertex>(pos3);
219 vertexmap2to3[*v2]=v3;
222 std::map<HepMC::GenParticle*,HepMC3::GenParticlePtr> particlemap2to3;
223 for (
auto p2=evt.particles_begin(); p2!=evt.particles_end(); p2++)
226 HepMC::FourVector mom2=(*p2)->momentum();
228 HepMC3::GenParticlePtr p3= std::make_shared<HepMC3::GenParticle>(mom3,(*p2)->pdg_id(),(*p2)->status());
230 p3->set_generated_mass((*p2)->generated_mass());
231 particlemap2to3[*p2]=p3;
233 auto v2production=(*p2)->production_vertex();
234 HepMC3::GenVertexPtr v3production;
237 auto v3=vertexmap2to3.find(v2production);
238 if (v3!=vertexmap2to3.end()) v3production=vertexmap2to3[v2production];
241 auto v2end=(*p2)->end_vertex();
242 HepMC3::GenVertexPtr v3end;
245 auto v3=vertexmap2to3.find(v2end);
246 if (v3!=vertexmap2to3.end()) v3end=vertexmap2to3[v2end];
249 if (v3production) v3production->add_particle_out(p3);
250 if (v3end) v3end->add_particle_in(p3);
252 std::shared_ptr<HepMC3::DoubleAttribute> p3_theta=std::make_shared<HepMC3::DoubleAttribute>((*p2)->polarization().theta());
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())
257 p3->add_attribute(
"theta",p3_theta);
258 p3->add_attribute(
"phi",p3_theta);
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);
266 p3->add_attribute(
"flows",flows);
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());
276 std::shared_ptr<HepMC3::GenPdfInfo> pdf3 = std::make_shared<HepMC3::GenPdfInfo>();
277 auto pdf2=evt.pdf_info();
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();
292 std::shared_ptr<HepMC3::GenHeavyIon> hi3 = std::make_shared<HepMC3::GenHeavyIon>();
293 auto hi2=evt.heavy_ion();
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();
315 std::stringstream ss;
316 evt.weights().print(ss);
317 std::string allweights=ss.str();
318 std::vector<std::string> wnames;
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;
329 n->
run_info()->set_weight_names(wnames);
330 n->
weights()=std::vector<double>(wnames.size(),1.0);
331 for (
auto wn: wnames)
333 n->
weight(wn)=evt.weights()[wn];
336 std::shared_ptr<HepMC3::DoubleAttribute> A_event_scale=std::make_shared<HepMC3::DoubleAttribute>(evt.event_scale());
338 std::shared_ptr<HepMC3::DoubleAttribute> A_alphaQED=std::make_shared<HepMC3::DoubleAttribute>(evt.alphaQED());
340 std::shared_ptr<HepMC3::DoubleAttribute> A_alphaQCD=std::make_shared<HepMC3::DoubleAttribute>(evt.alphaQCD());
342 std::shared_ptr<HepMC3::IntAttribute> A_signal_process_id=std::make_shared<HepMC3::IntAttribute>(evt.signal_process_id());
344 std::shared_ptr<HepMC3::IntAttribute> A_mpi=std::make_shared<HepMC3::IntAttribute>(evt.mpi());
346 std::shared_ptr<HepMC3::VectorLongIntAttribute> A_random_states= std::make_shared<HepMC3::VectorLongIntAttribute>( evt.random_states());
349 auto spv2=evt.signal_process_vertex();
350 if (vertexmap2to3.find(spv2)!=vertexmap2to3.end())
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);
const std::vector< ConstGenVertexPtr > & vertices() const
Get list of vertices (const)
void set_cross_section(GenCrossSectionPtr cs)
Set cross-section information.
void add_vertex(GenVertexPtr v)
Add vertex.
double t() const
Time component of position/displacement.
Definition of class GenParticle.
const std::vector< std::string > & weight_names() const
Definition of class GenVertex.
Attribute that holds a vector of integers of type int.
const Units::LengthUnit & length_unit() const
Get length unit.
double x() const
x-component of position/displacement
int event_number() const
Get event number.
void add_attribute(const std::string &name, const std::shared_ptr< Attribute > &att, const int &id=0)
Add event attribute to event.
void set_cross_section(const double &xs, const double &xs_err, const long &n_acc=-1, const long &n_att=-1)
Set all fields.
std::vector< int > value() const
get the value associated to this Attribute.
LengthUnit
Position units.
Stores additional information about PDFs.
std::vector< ConstGenParticlePtr > beams() const
Vector of beam particles.
MomentumUnit
Momentum units.
const Units::MomentumUnit & momentum_unit() const
Get momentum unit.
double e() const
Energy component of momentum.
Stores event-related information.
Attribute that holds a real number as a double.
void set_pdf_info(GenPdfInfoPtr pi)
Set PDF information.
double px() const
x-component of momentum
HepMC3::GenEvent * ConvertHepMCGenEvent_2to3(const HepMC::GenEvent &evt, std::shared_ptr< HepMC3::GenRunInfo > run)
std::shared_ptr< T > attribute(const std::string &name, const int &id=0) const
Get attribute of type T.
void set_run_info(std::shared_ptr< GenRunInfo > run)
Set the GenRunInfo object by smart pointer.
Stores additional information about Heavy Ion generator.
GenHeavyIon HeavyIon
Backward compatibility typedef.
void set_units(Units::MomentumUnit new_momentum_unit, Units::LengthUnit new_length_unit)
Change event units Converts event from current units to new ones.
HepMC::GenEvent * ConvertHepMCGenEvent_3to2(const HepMC3::GenEvent &evt)
Please note the HEPMC_HAS_CENTRALITY should be defined externaly.
double y() const
y-component of position/displacement
const std::vector< double > & weights() const
Get event weight values as a vector.
double pz() const
z-component of momentum
Definition of class GenEvent.
double py() const
y-component of momentum
double weight(const unsigned long &index=0) const
const std::vector< ConstGenParticlePtr > & particles() const
Get list of particles (const)
void set_heavy_ion(GenHeavyIonPtr hi)
Set heavy ion generator additional information.
int value() const
get the value associated to this Attribute.
std::shared_ptr< GenRunInfo > run_info() const
Get a pointer to the the GenRunInfo object.
Attribute that holds an Integer implemented as an int.
Feature< Feature_type > abs(const Feature< Feature_type > &input)
Obtain the absolute value of a Feature. This works as you'd expect. If foo is a valid Feature...
Stores additional information about cross-section.
void set_event_number(const int &num)
Set event number.