diff --git a/CMakeLists.txt b/CMakeLists.txt index 2a083ace097bfd46c265ff8eea54fe37b27be9bf..57fd2e167873a4969295eb7a251bae4d5ff344f3 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -19,7 +19,7 @@ else(CMAKE_BUILD_TYPE) endif(CMAKE_BUILD_TYPE) if (${CMAKE_VERSION} GREATER "3.8") -cmake_policy(SET CMP0069 NEW) + cmake_policy(SET CMP0069 NEW) endif() string(TOUPPER ${CMAKE_BUILD_TYPE} CMAKE_BUILD_TYPE) @@ -432,8 +432,8 @@ endif() # want to use the namespaced target names add_library(HepMC3::All ALIAS HepMC3All) if(HEPMC3_BUILD_STATIC_LIBS) -add_library(HepMC3::All_static ALIAS HepMC3All_static) -install(TARGETS HepMC3All HepMC3All_static EXPORT HepMC3Targets) + add_library(HepMC3::All_static ALIAS HepMC3All_static) + install(TARGETS HepMC3All HepMC3All_static EXPORT HepMC3Targets) else() install(TARGETS HepMC3All EXPORT HepMC3Targets) endif() diff --git a/include/HepMC3/GenCrossSection.h b/include/HepMC3/GenCrossSection.h index b1007ee21e9ae69aa8eae612de5a47009d7069f7..4bc21b09de833b3c88490f9f43f4188c1acbb94f 100644 --- a/include/HepMC3/GenCrossSection.h +++ b/include/HepMC3/GenCrossSection.h @@ -64,6 +64,18 @@ public: /** @brief Set all fields */ void set_cross_section(const double& xs, const double& xs_err,const long& n_acc = -1, const long& n_att = -1); + /** @brief Set all fields */ + void set_cross_section(const std::vector<double>& xs, const std::vector<double>& xs_err,const long& n_acc = -1, const long& n_att = -1); + + /** @brief Get the cross-sections + */ + const std::vector<double>& xsecs() const { return cross_sections; } + + /** @brief Get the cross-section errors + */ + const std::vector<double>& xsec_errs() const { return cross_section_errors; } + + /** @brief Set the number of accepted events */ void set_accepted_events(const long& n_acc ) { @@ -92,56 +104,70 @@ public: named \a wName. */ void set_xsec(const std::string& wName,const double& xs) { - set_xsec(windx(wName), xs); + int pos = windx(wName); + if ( pos < 0 ) throw std::runtime_error("GenCrossSection::set_xsec(const std::string&,const double&): no weight with given name in this run"); + set_xsec(pos, xs); } /** @brief Set the cross section corresponding to the weight with index \a indx. */ - void set_xsec(const int& indx, const double& xs) { - cross_sections[indx] = xs; + void set_xsec(const unsigned long& index, const double& xs) { + if ( index >= cross_sections.size() ) {throw std::runtime_error("GenCrossSection::set_xsec(const unsigned long&): index outside of range");} + cross_sections[index] = xs; } /** @brief Set the cross section error corresponding to the weight named \a wName. */ void set_xsec_err(const std::string& wName, const double& xs_err) { - set_xsec_err(windx(wName), xs_err); + int pos = windx(wName); + if ( pos < 0 ) throw std::runtime_error("GenCrossSection::set_xsec_err(const std::string&,const double&): no weight with given name in this run"); + set_xsec_err(pos, xs_err); } /** @brief Set the cross section error corresponding to the weight with index \a indx. */ - void set_xsec_err(const int& indx, const double& xs_err) { - cross_section_errors[indx] = xs_err; + void set_xsec_err(const unsigned long& index, const double& xs_err) { + if ( index >= cross_section_errors.size() ) {throw std::runtime_error("GenCrossSection::set_xsec_err(const unsigned long&): index outside of range");} + cross_section_errors[index] = xs_err; } /** @brief Get the cross section corresponding to the weight named \a wName. */ double xsec(const std::string& wName) const { - return xsec(windx(wName)); + int pos = windx(wName); + if ( pos < 0 ) throw std::runtime_error("GenCrossSection::xsec(const std::string&): no weight with given name in this run"); + return xsec(pos); } /** @brief Get the cross section corresponding to the weight with index \a indx. */ - double xsec(const int& indx = 0) const { - return cross_sections[indx]; + double xsec(const unsigned long& index = 0) const { + if ( index < cross_sections.size() ) { return cross_sections.at(index); } + else { throw std::runtime_error("GenCrossSection::xsec(const unsigned long&): index outside of range");} + return 0.0; } /** @brief Get the cross section error corresponding to the weight named \a wName. */ double xsec_err(const std::string& wName) const { - return xsec_err(windx(wName)); + int pos = windx(wName); + if ( pos < 0 ) throw std::runtime_error("GenCrossSection::xsec_err(const std::string&): no weight with given name in this run"); + return xsec_err(pos); } /** @brief Get the cross section error corresponding to the weight with index \a indx. */ - double xsec_err(const int& indx = 0) const { - return cross_section_errors[indx]; + double xsec_err(const unsigned long& index = 0) const { + if ( index < cross_section_errors.size() ) {return cross_section_errors.at(index);} + else { throw std::runtime_error("GenCrossSection::xsec_err(const unsigned long&): index outside of range");} + return 0.0; } bool operator==( const GenCrossSection& ) const; ///< Operator == diff --git a/python/CMakeLists.txt b/python/CMakeLists.txt index d9e3deb9d210f0203cd99383090cc903d1468b07..23b867ae76012d5266c54605c69eebf2a2153031 100644 --- a/python/CMakeLists.txt +++ b/python/CMakeLists.txt @@ -170,7 +170,7 @@ if(HEPMC3_ENABLE_PROTOBUFIO) configure_file("pyHepMC3.protobufIO.egg-info.in" "pyHepMC3.protobufIO.egg-info" IMMEDIATE @ONLY) endif() -set(autoBN ${CMAKE_CURRENT_SOURCE_DIR}/src/pyHepMC3.cpp ${CMAKE_CURRENT_SOURCE_DIR}/src/pyHepMC3_0.cpp ${CMAKE_CURRENT_SOURCE_DIR}/src/pyHepMC3_1.cpp +set(autoBN ${CMAKE_CURRENT_SOURCE_DIR}/src/pyHepMC3.cpp ${CMAKE_CURRENT_SOURCE_DIR}/src/pyHepMC3_0.cpp ${CMAKE_CURRENT_SOURCE_DIR}/src/pyHepMC3_1.cpp ${CMAKE_CURRENT_SOURCE_DIR}/src/pyHepMC3_2.cpp ${CMAKE_CURRENT_SOURCE_DIR}/src/pyHepMC3_3.cpp ${CMAKE_CURRENT_SOURCE_DIR}/src/pyHepMC3_4.cpp ${CMAKE_CURRENT_SOURCE_DIR}/src/pyHepMC3_5.cpp ${CMAKE_CURRENT_SOURCE_DIR}/src/pyHepMC3_6.cpp ${CMAKE_CURRENT_SOURCE_DIR}/src/pyHepMC3_7.cpp ${CMAKE_CURRENT_SOURCE_DIR}/src/pyHepMC3_8.cpp ${CMAKE_CURRENT_SOURCE_DIR}/src/pyHepMC3_9.cpp ${CMAKE_CURRENT_SOURCE_DIR}/src/pyHepMC3_10.cpp diff --git a/python/src/pyHepMC3_7.cpp b/python/src/pyHepMC3_7.cpp index 8549c306ad592219b63765d2b3114f8c417288dc..c66af54e4609056c9959205d88a822966f25ddab 100644 --- a/python/src/pyHepMC3_7.cpp +++ b/python/src/pyHepMC3_7.cpp @@ -100,20 +100,25 @@ void bind_pyHepMC3_7(std::function< pybind11::module &(std::string const &namesp cl.def("set_cross_section", [](HepMC3::GenCrossSection &o, const double & a0, const double & a1) -> void { return o.set_cross_section(a0, a1); }, "", pybind11::arg("xs"), pybind11::arg("xs_err")); cl.def("set_cross_section", [](HepMC3::GenCrossSection &o, const double & a0, const double & a1, const long & a2) -> void { return o.set_cross_section(a0, a1, a2); }, "", pybind11::arg("xs"), pybind11::arg("xs_err"), pybind11::arg("n_acc")); cl.def("set_cross_section", (void (HepMC3::GenCrossSection::*)(const double &, const double &, const long &, const long &)) &HepMC3::GenCrossSection::set_cross_section, "Set all fields \n\nC++: HepMC3::GenCrossSection::set_cross_section(const double &, const double &, const long &, const long &) --> void", pybind11::arg("xs"), pybind11::arg("xs_err"), pybind11::arg("n_acc"), pybind11::arg("n_att")); + cl.def("set_cross_section", [](HepMC3::GenCrossSection &o, const class std::vector<double> & a0, const class std::vector<double> & a1) -> void { return o.set_cross_section(a0, a1); }, "", pybind11::arg("xs"), pybind11::arg("xs_err")); + cl.def("set_cross_section", [](HepMC3::GenCrossSection &o, const class std::vector<double> & a0, const class std::vector<double> & a1, const long & a2) -> void { return o.set_cross_section(a0, a1, a2); }, "", pybind11::arg("xs"), pybind11::arg("xs_err"), pybind11::arg("n_acc")); + cl.def("set_cross_section", (void (HepMC3::GenCrossSection::*)(const class std::vector<double> &, const class std::vector<double> &, const long &, const long &)) &HepMC3::GenCrossSection::set_cross_section, "Set all fields \n\nC++: HepMC3::GenCrossSection::set_cross_section(const class std::vector<double> &, const class std::vector<double> &, const long &, const long &) --> void", pybind11::arg("xs"), pybind11::arg("xs_err"), pybind11::arg("n_acc"), pybind11::arg("n_att")); + cl.def("xsecs", (const class std::vector<double> & (HepMC3::GenCrossSection::*)() const) &HepMC3::GenCrossSection::xsecs, "Get the cross-sections\n\nC++: HepMC3::GenCrossSection::xsecs() const --> const class std::vector<double> &", pybind11::return_value_policy::automatic); + cl.def("xsec_errs", (const class std::vector<double> & (HepMC3::GenCrossSection::*)() const) &HepMC3::GenCrossSection::xsec_errs, "Get the cross-section errors\n\nC++: HepMC3::GenCrossSection::xsec_errs() const --> const class std::vector<double> &", pybind11::return_value_policy::automatic); cl.def("set_accepted_events", (void (HepMC3::GenCrossSection::*)(const long &)) &HepMC3::GenCrossSection::set_accepted_events, "Set the number of accepted events\n\nC++: HepMC3::GenCrossSection::set_accepted_events(const long &) --> void", pybind11::arg("n_acc")); cl.def("set_attempted_events", (void (HepMC3::GenCrossSection::*)(const long &)) &HepMC3::GenCrossSection::set_attempted_events, "Set the number of attempted events\n\nC++: HepMC3::GenCrossSection::set_attempted_events(const long &) --> void", pybind11::arg("n_att")); cl.def("get_accepted_events", (long (HepMC3::GenCrossSection::*)() const) &HepMC3::GenCrossSection::get_accepted_events, "Get the number of accepted events\n\nC++: HepMC3::GenCrossSection::get_accepted_events() const --> long"); cl.def("get_attempted_events", (long (HepMC3::GenCrossSection::*)() const) &HepMC3::GenCrossSection::get_attempted_events, "Get the number of attempted events\n\nC++: HepMC3::GenCrossSection::get_attempted_events() const --> long"); cl.def("set_xsec", (void (HepMC3::GenCrossSection::*)(const std::string &, const double &)) &HepMC3::GenCrossSection::set_xsec, "Set the cross section corresponding to the weight\n named \n \n\nC++: HepMC3::GenCrossSection::set_xsec(const std::string &, const double &) --> void", pybind11::arg("wName"), pybind11::arg("xs")); - cl.def("set_xsec", (void (HepMC3::GenCrossSection::*)(const int &, const double &)) &HepMC3::GenCrossSection::set_xsec, "Set the cross section corresponding to the weight with\n index \n \n\nC++: HepMC3::GenCrossSection::set_xsec(const int &, const double &) --> void", pybind11::arg("indx"), pybind11::arg("xs")); + cl.def("set_xsec", (void (HepMC3::GenCrossSection::*)(const unsigned long &, const double &)) &HepMC3::GenCrossSection::set_xsec, "Set the cross section corresponding to the weight with\n index \n \n\nC++: HepMC3::GenCrossSection::set_xsec(const unsigned long &, const double &) --> void", pybind11::arg("index"), pybind11::arg("xs")); cl.def("set_xsec_err", (void (HepMC3::GenCrossSection::*)(const std::string &, const double &)) &HepMC3::GenCrossSection::set_xsec_err, "Set the cross section error corresponding to the weight\n named \n \n\nC++: HepMC3::GenCrossSection::set_xsec_err(const std::string &, const double &) --> void", pybind11::arg("wName"), pybind11::arg("xs_err")); - cl.def("set_xsec_err", (void (HepMC3::GenCrossSection::*)(const int &, const double &)) &HepMC3::GenCrossSection::set_xsec_err, "Set the cross section error corresponding to the weight\n with index \n \n\nC++: HepMC3::GenCrossSection::set_xsec_err(const int &, const double &) --> void", pybind11::arg("indx"), pybind11::arg("xs_err")); + cl.def("set_xsec_err", (void (HepMC3::GenCrossSection::*)(const unsigned long &, const double &)) &HepMC3::GenCrossSection::set_xsec_err, "Set the cross section error corresponding to the weight\n with index \n \n\nC++: HepMC3::GenCrossSection::set_xsec_err(const unsigned long &, const double &) --> void", pybind11::arg("index"), pybind11::arg("xs_err")); cl.def("xsec", (double (HepMC3::GenCrossSection::*)(const std::string &) const) &HepMC3::GenCrossSection::xsec, "Get the cross section corresponding to the weight named\n \n \n\nC++: HepMC3::GenCrossSection::xsec(const std::string &) const --> double", pybind11::arg("wName")); cl.def("xsec", [](HepMC3::GenCrossSection const &o) -> double { return o.xsec(); }, ""); - cl.def("xsec", (double (HepMC3::GenCrossSection::*)(const int &) const) &HepMC3::GenCrossSection::xsec, "Get the cross section corresponding to the weight with index\n \n \n\nC++: HepMC3::GenCrossSection::xsec(const int &) const --> double", pybind11::arg("indx")); + cl.def("xsec", (double (HepMC3::GenCrossSection::*)(const unsigned long &) const) &HepMC3::GenCrossSection::xsec, "Get the cross section corresponding to the weight with index\n \n \n\nC++: HepMC3::GenCrossSection::xsec(const unsigned long &) const --> double", pybind11::arg("index")); cl.def("xsec_err", (double (HepMC3::GenCrossSection::*)(const std::string &) const) &HepMC3::GenCrossSection::xsec_err, "Get the cross section error corresponding to the weight\n named \n \n\nC++: HepMC3::GenCrossSection::xsec_err(const std::string &) const --> double", pybind11::arg("wName")); cl.def("xsec_err", [](HepMC3::GenCrossSection const &o) -> double { return o.xsec_err(); }, ""); - cl.def("xsec_err", (double (HepMC3::GenCrossSection::*)(const int &) const) &HepMC3::GenCrossSection::xsec_err, "Get the cross section error corresponding to the weight\n with index \n \n\nC++: HepMC3::GenCrossSection::xsec_err(const int &) const --> double", pybind11::arg("indx")); + cl.def("xsec_err", (double (HepMC3::GenCrossSection::*)(const unsigned long &) const) &HepMC3::GenCrossSection::xsec_err, "Get the cross section error corresponding to the weight\n with index \n \n\nC++: HepMC3::GenCrossSection::xsec_err(const unsigned long &) const --> double", pybind11::arg("index")); cl.def("__eq__", (bool (HepMC3::GenCrossSection::*)(const class HepMC3::GenCrossSection &) const) &HepMC3::GenCrossSection::operator==, "C++: HepMC3::GenCrossSection::operator==(const class HepMC3::GenCrossSection &) const --> bool", pybind11::arg("")); cl.def("__ne__", (bool (HepMC3::GenCrossSection::*)(const class HepMC3::GenCrossSection &) const) &HepMC3::GenCrossSection::operator!=, "C++: HepMC3::GenCrossSection::operator!=(const class HepMC3::GenCrossSection &) const --> bool", pybind11::arg("")); cl.def("is_valid", (bool (HepMC3::GenCrossSection::*)() const) &HepMC3::GenCrossSection::is_valid, "C++: HepMC3::GenCrossSection::is_valid() const --> bool"); diff --git a/src/GenCrossSection.cc b/src/GenCrossSection.cc index 0cb8c56eb3ac62973f07be542a97683c5ad4028c..77f462acf3dc2d4755eaa4921c39996910f37e14 100644 --- a/src/GenCrossSection.cc +++ b/src/GenCrossSection.cc @@ -29,13 +29,20 @@ void GenCrossSection::set_cross_section(const double& xs, const double& xs_err, double cross_section_error = xs_err; accepted_events = n_acc; attempted_events = n_att; - size_t N = 1; - if ( event() ) N = std::max(event()->weights().size(), N); + size_t N = std::max( event() ? event()->weights().size() : 0, size_t{1}); cross_sections = std::vector<double>(N, cross_section); cross_section_errors = std::vector<double>(N, cross_section_error); } +void GenCrossSection::set_cross_section(const std::vector<double>& xs, const std::vector<double>& xs_err, const long& n_acc, const long& n_att) { + cross_sections = xs; + cross_section_errors = xs_err; + accepted_events = n_acc; + attempted_events = n_att; +} + + bool GenCrossSection::from_string(const std::string &att) { const char *cursor = att.data(); cross_sections.clear();