Page MenuHomeHEPForge
Files F8723710

No OneTemporary
Actions

View Options

diff --git a/include/RHEJ/HepMCInterface.hh b/include/RHEJ/HepMCInterface.hh
new file mode 100644
index 0000000..00fb8aa
--- /dev/null
+++ b/include/RHEJ/HepMCInterface.hh
@@ -0,0 +1,44 @@
+/** \file HepMCInterface.hh
+ * \brief Converts HEJ events into HepMC::GenEvents
+ */
+
+#ifdef RHEJ_BUILD_WITH_HepMC_VERSION
+#pragma once
+
+#include <cstddef>
+
+#if RHEJ_BUILD_WITH_HepMC_VERSION >= 3
+#include "HepMC/GenRunInfo.h"
+#endif
+
+namespace LHEF{
+ struct HEPRUP;
+}
+namespace HepMC{
+ class GenEvent;
+ class GenCrossSection;
+}
+
+namespace RHEJ{
+ class Event;
+
+ class HepMCInterface{
+ public:
+ HepMCInterface(LHEF::HEPRUP &);
+ HepMC::GenEvent operator()(Event const & event, size_t const weight_index = 0);
+ private:
+ size_t event_count;
+ double tot_weight;
+ double tot_weight2;
+ public:
+ void add_event(double const wt);
+ #if RHEJ_BUILD_WITH_HepMC_VERSION >= 3
+ HepMC::shared_ptr<HepMC::GenRunInfo> runinfo;
+ HepMC::shared_ptr<HepMC::GenCrossSection> cross_section();
+ #else
+ HepMC::GenCrossSection cross_section();
+ #endif
+ };
+}
+
+#endif
diff --git a/src/HepMCInteface.cc b/src/HepMCInteface.cc
new file mode 100644
index 0000000..95e9a79
--- /dev/null
+++ b/src/HepMCInteface.cc
@@ -0,0 +1,171 @@
+#include "RHEJ/HepMCInterface.hh"
+
+#ifdef RHEJ_BUILD_WITH_HepMC_VERSION
+
+#include "RHEJ/Event.hh"
+#include "RHEJ/Version.hh"
+
+#include "HepMC/GenEvent.h"
+#include "HepMC/GenVertex.h"
+#include "HepMC/GenParticle.h"
+#include "HepMC/GenCrossSection.h"
+
+#if RHEJ_BUILD_WITH_HepMC_VERSION >= 3
+#include "HepMC/LHEFAttributes.h"
+#endif
+
+#include "LHEF/LHEF.h"
+
+namespace RHEJ{
+
+ namespace {
+ HepMC::FourVector to_FourVector(Sparticle const & sp){
+ return {sp.px(), sp.py(), sp.pz(), sp.E()};
+ }
+ #if RHEJ_BUILD_WITH_HepMC_VERSION >= 3
+ void add_generator_tag(LHEF::HEPRUP & heprup){
+ heprup.generators.emplace_back(LHEF::XMLTag{});
+ heprup.generators.back().name = Version::package_name();
+ heprup.generators.back().version = Version::String();
+ }
+ void reset_weight_info(LHEF::HEPRUP & heprup){
+ heprup.IDWTUP = -4;
+ // use placeholders for unknown init block values
+ // we can overwrite them after processing all events
+ heprup.XSECUP = {0.};
+ heprup.XERRUP = {0.};
+ heprup.XMAXUP = {0.};
+ }
+
+ HepMC::shared_ptr<HepMC::GenRunInfo> & add_HEPRUP(
+ HepMC::shared_ptr<HepMC::GenRunInfo> & runinfo, LHEF::HEPRUP heprup
+ ){
+ //@TODO something in here gives a stackfault
+ add_generator_tag(heprup);
+ reset_weight_info(heprup);
+ auto hepr = HepMC::make_shared<HepMC::HEPRUPAttribute>();
+ hepr->heprup = std::move(heprup);
+ runinfo->add_attribute(std::string("HEPRUP"), hepr);
+
+ for (int i = 0, N = hepr->heprup.generators.size(); i < N; ++i ){
+ HepMC::GenRunInfo::ToolInfo tool;
+ tool.name = hepr->heprup.generators[i].name;
+ tool.version = hepr->heprup.generators[i].version;
+ tool.description = hepr->heprup.generators[i].contents;
+ runinfo->tools().push_back(tool);
+ }
+ return runinfo;
+ }
+ #endif // HepMC 3
+ } // namespace anonymous
+
+#if RHEJ_BUILD_WITH_HepMC_VERSION >= 3
+ HepMCInterface::HepMCInterface(LHEF::HEPRUP & heprup):
+ event_count(0.), tot_weight(0.), tot_weight2(0.),
+ runinfo(add_HEPRUP(runinfo, std::move(heprup)))
+ {}
+#else
+ HepMCInterface::HepMCInterface(LHEF::HEPRUP &):
+ event_count(0.), tot_weight(0.), tot_weight2(0.)
+ {}
+#endif
+
+ void HepMCInterface::add_event(double const wt){
+ ++event_count;
+ tot_weight += wt;
+ tot_weight2 += wt * wt;
+ }
+
+#if RHEJ_BUILD_WITH_HepMC_VERSION >= 3
+ HepMC::shared_ptr<HepMC::GenCrossSection> HepMCInterface::cross_section(){
+ auto xs = HepMC::make_shared<HepMC::GenCrossSection>();
+ xs->set_cross_section(tot_weight, sqrt(tot_weight2));
+ return xs;
+ }
+#else
+ HepMC::GenCrossSection HepMCInterface::cross_section(){
+ HepMC::GenCrossSection xs;
+ xs.set_cross_section(tot_weight, sqrt(tot_weight2));
+ return xs;
+ }
+#endif
+
+#if RHEJ_BUILD_WITH_HepMC_VERSION >= 3 // TODO better way?
+ #define array_type(x) HepMC::make_shared<x>
+#else // HepMC 2
+ #define array_type(x) new x
+#endif
+
+ HepMC::GenEvent HepMCInterface::operator()(Event const & event, size_t const weight_index){
+ HepMC::GenEvent out_ev{HepMC::Units::GEV, HepMC::Units::MM};
+ auto vx = array_type(HepMC::GenVertex)();
+ for(auto const & in: event.incoming()){
+ vx->add_particle_in(
+ array_type(HepMC::GenParticle)(
+ to_FourVector(in), static_cast<int>(in.type), -1
+ )
+ );
+ }
+ for(int i=0; i< (int) event.outgoing().size(); ++i){
+ auto const & out = event.outgoing()[i];
+ auto particle = array_type(HepMC::GenParticle)(
+ to_FourVector(out), static_cast<int>(out.type), +1
+ );
+ for(auto const & decay: event.decays()){
+ if( i == decay.first ){
+ particle->set_status(+3);
+ auto vx_decay= array_type(HepMC::GenVertex)();
+ vx_decay->add_particle_in(particle);
+ for( auto const & out: decay.second){
+ vx_decay->add_particle_out(
+ array_type(HepMC::GenParticle)(
+ to_FourVector(out), static_cast<int>(out.type), +1
+ )
+ );
+ }
+ out_ev.add_vertex(vx_decay);
+ break;
+ }
+ }
+ vx->add_particle_out(particle);
+ }
+ out_ev.add_vertex(vx);
+ /// weights
+
+ std::cout << event.variations().size() << std::endl;
+ EventParameters central_parameters;
+ if(event.variations().size() == 0){
+ central_parameters = event.central();
+ }
+ else{
+ central_parameters = event.variations(weight_index);
+ }
+ out_ev.weights().push_back(central_parameters.weight);
+ for(auto const & var: event.variations()){
+ out_ev.weights().push_back(var.weight);
+ }
+ /// @TODO add name list
+
+ #if RHEJ_BUILD_WITH_HepMC_VERSION >= 3
+ std::cout << "Hep MC 3" << std::endl;
+ #else // HepMC 2
+ std::cout << "Hep MC 2" << std::endl;
+ out_ev.set_signal_process_id(event.type()+1); // "+1": conistent with lhe
+ out_ev.set_event_scale(central_parameters.mur); // event scale
+ #endif
+ /// cross section
+ add_event(central_parameters.weight);
+
+ // out_ev.set_event_number(genInfo_.event_count); /// event number
+
+ out_ev.set_cross_section( cross_section() );
+
+ /// @TODO add alphaQCD (need function) and alphaQED
+ /// @TODO output pdf (currently not avaiable from event alone)
+ return out_ev;
+ }
+ #undef array_type
+
+
+}
+#endif // HepMC
diff --git a/src/HepMCWriter.cc b/src/HepMCWriter.cc
index 9c4fd34..c2f003d 100644
--- a/src/HepMCWriter.cc
+++ b/src/HepMCWriter.cc
@@ -1,251 +1,118 @@
#include "RHEJ/HepMCWriter.hh"
#include <cassert>
#ifdef RHEJ_BUILD_WITH_HepMC_VERSION
#if RHEJ_BUILD_WITH_HepMC_VERSION >= 3
#include "HepMC/WriterAscii.h"
-#include "HepMC/LHEFAttributes.h"
#else
#include "HepMC/IO_GenEvent.h"
#endif
#include "HepMC/GenVertex.h"
#include "HepMC/GenParticle.h"
#include "RHEJ/utility.hh"
#include "RHEJ/Event.hh"
#include "RHEJ/exceptions.hh"
-#include "RHEJ/Version.hh"
+#include "RHEJ/HepMCInterface.hh"
namespace RHEJ{
#if RHEJ_BUILD_WITH_HepMC_VERSION >= 3
- namespace {
- void add_generator_tag(LHEF::HEPRUP & heprup){
- heprup.generators.emplace_back(LHEF::XMLTag{});
- heprup.generators.back().name = Version::package_name();
- heprup.generators.back().version = Version::String();
- }
- void reset_weight_info(LHEF::HEPRUP & heprup){
- heprup.IDWTUP = -4;
- // use placeholders for unknown init block values
- // we can overwrite them after processing all events
- heprup.XSECUP = {0.};
- heprup.XERRUP = {0.};
- heprup.XMAXUP = {0.};
- }
- HepMC::FourVector to_FourVector(Sparticle const & sp){
- return {sp.px(), sp.py(), sp.pz(), sp.E()};
- }
-
- HepMC::shared_ptr<HepMC::GenRunInfo> & add_HEPRUP(
- HepMC::shared_ptr<HepMC::GenRunInfo> & runinfo,
- LHEF::HEPRUP heprup
- ){
- add_generator_tag(heprup);
- reset_weight_info(heprup);
-
- auto hepr = HepMC::make_shared<HepMC::HEPRUPAttribute>();
- hepr->heprup = std::move(heprup);
- runinfo->add_attribute("HEPRUP", hepr);
-
- for (int i = 0, N = hepr->heprup.generators.size(); i < N; ++i ){
- HepMC::GenRunInfo::ToolInfo tool;
- tool.name = hepr->heprup.generators[i].name;
- tool.version = hepr->heprup.generators[i].version;
- tool.description = hepr->heprup.generators[i].contents;
- runinfo->tools().push_back(tool);
- }
- return runinfo;
- }
- } // namespace anonymous
struct HepMCWriter::HepMCWriterImpl{
-
- HepMC::shared_ptr<HepMC::GenRunInfo> runinfo_;
+ HepMCInterface hepmc_;
HepMC::WriterAscii writer_;
HepMCWriterImpl(
std::string const & file, LHEF::HEPRUP && heprup
):
- runinfo_{HepMC::make_shared<HepMC::GenRunInfo>()},
- writer_{file, add_HEPRUP(runinfo_, std::move(heprup))}
+ hepmc_(heprup),
+ writer_{file, hepmc_.runinfo}
{}
HepMCWriterImpl & operator=(HepMCWriterImpl const & other) = delete;
HepMCWriterImpl(HepMCWriterImpl const & other) = delete;
HepMCWriterImpl & operator=(HepMCWriterImpl && other) = delete;
HepMCWriterImpl(HepMCWriterImpl && other) = delete;
~HepMCWriterImpl(){
writer_.close();
}
void write(Event const & ev){
- if(!ev.decays().empty()){
- throw not_implemented{"HepMC output for events with decays"};
- }
- auto vx = HepMC::make_shared<HepMC::GenVertex>();
- for(auto const & in: ev.incoming()){
- vx->add_particle_in(
- HepMC::make_shared<HepMC::GenParticle>(
- to_FourVector(in), static_cast<int>(in.type), -1
- )
- );
- }
- for(auto const & out: ev.outgoing()){
- vx->add_particle_out(
- HepMC::make_shared<HepMC::GenParticle>(
- to_FourVector(out), static_cast<int>(out.type), +1
- )
- );
- }
- HepMC::GenEvent out_ev{runinfo_, HepMC::Units::GEV, HepMC::Units::MM};
- out_ev.add_vertex(vx);
- out_ev.weights().reserve(ev.variations().size() + 1u);
- out_ev.weights().emplace_back(ev.central().weight);
- for(auto const & var: ev.variations()){
- out_ev.weights().emplace_back(var.weight);
- }
+ auto out_ev = hepmc_(ev);
writer_.write_event(out_ev);
}
};
void HepMCWriter::HepMCWriterImplDeleter::operator()(HepMCWriterImpl* p) {
delete p;
}
#else // HepMC 2
- namespace{
- HepMC::FourVector to_FourVector(Sparticle const & sp){
- return {sp.px(), sp.py(), sp.pz(), sp.E()};
- }
- }
-
struct HepMCWriter::HepMCWriterImpl{
- struct HepMCInfo{
- size_t event_count = 0.;
- double tot_weight = 0.;
- double tot_weight2 = 0.;
-
- void add_event(Event const & ev){
- double wt = ev.central().weight;
- tot_weight += wt;
- tot_weight2 += wt * wt;
- ++event_count;
- }
- };
-
HepMC::IO_GenEvent writer_;
- HepMCInfo genInfo_;
+ HepMCInterface hepmc_;
HepMCWriterImpl(
- std::string const & file, LHEF::HEPRUP &&
+ std::string const & file, LHEF::HEPRUP && heprup
):
- writer_{file}
+ writer_{file},
+ hepmc_(heprup)
{}
HepMCWriterImpl & operator=(HepMCWriterImpl const & other) = delete;
HepMCWriterImpl(HepMCWriterImpl const & other) = delete;
HepMCWriterImpl & operator=(HepMCWriterImpl && other) = delete;
HepMCWriterImpl(HepMCWriterImpl && other) = delete;
void write(Event const & ev){
- auto vx = new HepMC::GenVertex();
- for(auto const & in: ev.incoming()){
- vx->add_particle_in(
- new HepMC::GenParticle{
- to_FourVector(in), static_cast<int>(in.type), -1
- }
- );
- }
- for(int i=0; i< (int) ev.outgoing().size(); ++i){
- auto const & out = ev.outgoing()[i];
- bool decays = false;
- for(auto const & decay: ev.decays()){
- if( i == decay.first ){ // if particle decays replace it with decay products
- for( auto const & out: decay.second){
- vx->add_particle_out(
- new HepMC::GenParticle{
- to_FourVector(out), static_cast<int>(out.type), +1
- }
- );
- }
- decays = true;
- break;
- }
- }
- if(!decays)
- vx->add_particle_out(
- new HepMC::GenParticle{
- to_FourVector(out), static_cast<int>(out.type), +1
- }
- );
- }
- HepMC::GenEvent out_ev{HepMC::Units::GEV, HepMC::Units::MM};
- out_ev.add_vertex(vx);
- /// weights
- out_ev.weights().push_back(ev.central().weight);
- for(auto const & var: ev.variations()){
- out_ev.weights().push_back(var.weight);
- }
- /// @TODO add name list
- /// general informations
- genInfo_.add_event(ev);
- out_ev.set_signal_process_id(ev.type()+1); /// "+1": conistent with lhe
- out_ev.set_event_scale(ev.central().mur); /// event scale
- out_ev.set_event_number(genInfo_.event_count); /// event number
- /// cross section
- HepMC::GenCrossSection t_xs;
- t_xs.set_cross_section( genInfo_.tot_weight , sqrt(genInfo_.tot_weight2) );
- out_ev.set_cross_section( t_xs );
-
- /// @TODO add alphaQCD (need function) and alphaQED
- /// @TODO output pdf (currently not avaiable from event alone)
+ auto out_ev = hepmc_(ev);
writer_.write_event(&out_ev);
}
};
void HepMCWriter::HepMCWriterImplDeleter::operator()(HepMCWriterImpl* p) {
delete p;
}
#endif // HepMC 2
HepMCWriter::HepMCWriter(std::string const & file, LHEF::HEPRUP heprup):
impl_{std::unique_ptr<HepMCWriterImpl, HepMCWriterImplDeleter>{
new HepMCWriterImpl(file, std::move(heprup))
}}
{}
void HepMCWriter::write(Event const & ev){
impl_->write(ev);
}
} // namespace RHEJ
#else // no HepMC
namespace RHEJ{
class HepMCWriter::HepMCWriterImpl{};
HepMCWriter::HepMCWriter(std::string const &, LHEF::HEPRUP){
throw std::invalid_argument(
"Failed to create HepMC writer: "
"Reversed HEJ was built without HepMC support"
);
}
void HepMCWriter::write(Event const &){
assert(false);
}
void HepMCWriter::HepMCWriterImplDeleter::operator()(HepMCWriterImpl* p) {
delete p;
}
}
#endif

File Metadata

Mime Type
text/x-diff
Expires
Mon, Jan 20, 9:12 PM (1 d, 1 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
4242419
Default Alt Text
(14 KB)

Event Timeline