Page MenuHomeHEPForge
Files F10664257

No OneTemporary
Actions

View Options

diff --git a/FixedOrderGen/src/config.cc b/FixedOrderGen/src/config.cc
index fe4b9dc..38a4c88 100644
--- a/FixedOrderGen/src/config.cc
+++ b/FixedOrderGen/src/config.cc
@@ -1,305 +1,302 @@
#include "config.hh"
#include <cctype>
#include "HEJ/config.hh"
#include "HEJ/YAMLreader.hh"
namespace HEJFOG{
using HEJ::set_from_yaml;
using HEJ::set_from_yaml_if_defined;
namespace{
//! Get YAML tree of supported options
/**
* The configuration file is checked against this tree of options
* in assert_all_options_known.
*/
YAML::Node const & get_supported_options(){
const static YAML::Node supported = [](){
YAML::Node supported;
static const auto opts = {
"process", "events", "unordered fraction", "scales", "scale factors",
"max scale ratio", "pdf", "event output", "analysis", "import scales"
};
// add subnodes to "supported" - the assigned value is irrelevant
for(auto && opt: opts) supported[opt] = "";
for(auto && jet_opt: {"min pt", "peak pt", "algorithm", "R", "max rapidity"}){
supported["jets"][jet_opt] = "";
}
for(auto && particle_type: {"Higgs", "Wp", "Wm", "Z"}){
for(auto && particle_opt: {"mass", "width"}){
supported["Particle properties"][particle_type][particle_opt] = "";
}
supported["Particle properties"][particle_type]["decays"]["into"] = "";
supported["Particle properties"][particle_type]["decays"]["branching ratio"] = "";
}
for(auto && opt: {"mt", "use impact factors", "include bottom", "mb"}){
supported["Higgs coupling"][opt] = "";
}
for(auto && beam_opt: {"energy", "particles"}){
supported["beam"][beam_opt] = "";
}
for(auto && unweight_opt: {"sample size", "max deviation"}){
supported["unweight"][unweight_opt] = "";
}
for(auto && opt: {"name", "seed"}){
supported["random generator"][opt] = "";
}
return supported;
}();
return supported;
}
JetParameters get_jet_parameters(
YAML::Node const & node, std::string const & entry
){
const auto p = HEJ::get_jet_parameters(node, entry);
JetParameters result;
result.def = p.def;
result.min_pt = p.min_pt;
set_from_yaml(result.max_y, node, entry, "max rapidity");
set_from_yaml_if_defined(result.peak_pt, node, entry, "peak pt");
return result;
}
Beam get_Beam(
YAML::Node const & node, std::string const & entry
){
Beam beam;
std::vector<HEJ::ParticleID> particles;
set_from_yaml(beam.energy, node, entry, "energy");
set_from_yaml_if_defined(particles, node, entry, "particles");
if(! particles.empty()){
for(HEJ::ParticleID particle: particles){
if(particle != HEJ::pid::p && particle != HEJ::pid::p_bar){
throw std::invalid_argument{
"Unsupported value in option " + entry + ": particles:"
" only proton ('p') and antiproton ('p_bar') beams are supported"
};
}
}
if(particles.size() != 2){
throw std::invalid_argument{"Not exactly two beam particles"};
}
beam.particles.front() = particles.front();
beam.particles.back() = particles.back();
}
return beam;
}
std::vector<std::string> split(
std::string const & str, std::string const & delims
){
std::vector<std::string> result;
for(size_t begin, end = 0; end != str.npos;){
begin = str.find_first_not_of(delims, end);
if(begin == str.npos) break;
end = str.find_first_of(delims, begin + 1);
result.emplace_back(str.substr(begin, end - begin));
}
return result;
}
std::invalid_argument invalid_incoming(std::string const & what){
return std::invalid_argument{
"Incoming particle type " + what + " not supported,"
" incoming particles have to be 'p', 'p_bar' or partons"
};
}
std::invalid_argument invalid_outgoing(std::string const & what){
return std::invalid_argument{
"Outgoing particle type " + what + " not supported,"
" outgoing particles have to be 'j', 'photon', 'W+', 'W-', 'Z', 'H'"
};
}
Process get_process(
YAML::Node const & node, std::string const & entry
){
Process result;
std::string process_string;
set_from_yaml(process_string, node, entry);
assert(! process_string.empty());
const auto particles = split(process_string, " \n\t\v=>");
if(particles.size() < 3){
throw std::invalid_argument{
"Bad format in option process: '" + process_string
+ "', expected format is 'in1 in2 => out1 ...'"
};
}
result.incoming.front() = HEJ::to_ParticleID(particles[0]);
result.incoming.back() = HEJ::to_ParticleID(particles[1]);
for(size_t i = 0; i < result.incoming.size(); ++i){
const HEJ::ParticleID in = result.incoming[i];
if(
in != HEJ::pid::proton && in != HEJ::pid::p_bar
&& !HEJ::is_parton(in)
){
throw invalid_incoming(particles[i]);
}
}
result.njets = 0;
for(size_t i = result.incoming.size(); i < particles.size(); ++i){
assert(! particles[i].empty());
if(particles[i] == "j") ++result.njets;
- else if(std::isdigit(particles[i].front())){
- if(particles[i].back() != 'j'){
- throw invalid_outgoing(particles[i]);
- }
+ else if(std::isdigit(particles[i].front())
+ && particles[i].back() == 'j')
result.njets += std::stoi(particles[i]);
- }
else{
const auto pid = HEJ::to_ParticleID(particles[i]);
if(!HEJ::is_AWZH_boson(pid)){
throw invalid_outgoing(particles[i]);
}
if(result.boson){
throw std::invalid_argument{
"More than one outgoing boson is not supported"
};
}
result.boson = pid;
}
}
if(result.njets < 2){
throw std::invalid_argument{
"Process has to include at least two jets ('j')"
};
}
return result;
}
Decay get_decay(YAML::Node const & node){
Decay decay;
set_from_yaml(decay.products, node, "into");
set_from_yaml(decay.branching_ratio, node, "branching ratio");
return decay;
}
std::vector<Decay> get_decays(YAML::Node const & node){
using YAML::NodeType;
if(!node) return {};
switch(node.Type()){
case NodeType::Null:
case NodeType::Undefined:
return {};
case NodeType::Scalar:
throw HEJ::invalid_type{"value is not a list of decays"};
case NodeType::Map:
return {get_decay(node)};
case NodeType::Sequence:
std::vector<Decay> result;
for(auto && decay_str: node){
result.emplace_back();
set_from_yaml(result.back().products, decay_str, "into");
set_from_yaml(result.back().branching_ratio, decay_str, "branching ratio");
}
return result;
}
throw std::logic_error{"unreachable"};
}
ParticleProperties get_particle_properties(
YAML::Node const & node, std::string const & entry
){
ParticleProperties result;
set_from_yaml(result.mass, node, entry, "mass");
set_from_yaml(result.width, node, entry, "width");
try{
result.decays = get_decays(node[entry]["decays"]);
}
catch(HEJ::missing_option const & ex){
throw HEJ::missing_option{entry + ": decays: " + ex.what()};
}
catch(HEJ::invalid_type const & ex){
throw HEJ::invalid_type{entry + ": decays: " + ex.what()};
}
return result;
}
ParticlesPropMap get_all_particles_properties(YAML::Node const & node){
ParticlesPropMap result;
using namespace HEJ;
// @TODO allow more synonyms
if(node["Higgs"])
result[pid::Higgs] = get_particle_properties(node,"Higgs");
if(node["Wp"])
result[pid::Wp] = get_particle_properties(node,"Wp");
if(node["Wm"])
result[pid::Wm] = get_particle_properties(node,"Wm");
if(node["Z"])
result[pid::Z] = get_particle_properties(node,"Z");
return result;
}
UnweightSettings get_unweight(
YAML::Node const & node, std::string const & entry
){
UnweightSettings result;
set_from_yaml(result.sample_size, node, entry, "sample size");
if(result.sample_size <= 0){
throw std::invalid_argument{
"negative sample size " + std::to_string(result.sample_size)
};
}
set_from_yaml(result.max_dev, node, entry, "max deviation");
return result;
}
Config to_Config(YAML::Node const & yaml){
try{
HEJ::assert_all_options_known(yaml, get_supported_options());
}
catch(HEJ::unknown_option const & ex){
throw HEJ::unknown_option{std::string{"Unknown option '"} + ex.what() + "'"};
}
Config config;
config.process = get_process(yaml, "process");
set_from_yaml(config.events, yaml, "events");
config.jets = get_jet_parameters(yaml, "jets");
config.beam = get_Beam(yaml, "beam");
for(size_t i = 0; i < config.process.incoming.size(); ++i){
const auto & in = config.process.incoming[i];
using namespace HEJ::pid;
if( (in == p || in == p_bar) && in != config.beam.particles[i]){
throw std::invalid_argument{
"Particle type of beam " + std::to_string(i+1) + " incompatible"
+ " with type of incoming particle " + std::to_string(i+1)
};
}
}
set_from_yaml(config.pdf_id, yaml, "pdf");
set_from_yaml(config.unordered_fraction, yaml, "unordered fraction");
if(config.unordered_fraction < 0 || config.unordered_fraction > 1){
throw std::invalid_argument{
"unordered fraction has to be between 0 and 1"
};
}
config.particles_properties = get_all_particles_properties(
yaml["Particle properties"]);
set_from_yaml_if_defined(config.analysis_parameters, yaml, "analysis");
config.scales = HEJ::to_ScaleConfig(yaml);
set_from_yaml_if_defined(config.output, yaml, "event output");
config.rng = HEJ::to_RNGConfig(yaml, "random generator");
config.Higgs_coupling = HEJ::get_Higgs_coupling(yaml, "Higgs coupling");
if(yaml["unweight"]) config.unweight = get_unweight(yaml, "unweight");
return config;
}
} // namespace anonymous
Config load_config(std::string const & config_file){
try{
return to_Config(YAML::LoadFile(config_file));
}
catch(...){
std::cerr << "Error reading " << config_file << ":\n ";
throw;
}
}
}
diff --git a/include/HEJ/PDG_codes.hh b/include/HEJ/PDG_codes.hh
index 8b93ca3..01dc521 100644
--- a/include/HEJ/PDG_codes.hh
+++ b/include/HEJ/PDG_codes.hh
@@ -1,92 +1,97 @@
/** \file PDG_codes.hh
* \brief Contains the Particle IDs of all relevant SM particles.
*
* Large enumeration included which has multiple entries for potential
* alternative names of different particles. There are also functions
* which can be used to determine if a particle is a parton or if
* it is a non-gluon boson.
*/
#pragma once
+#include <string>
+
namespace HEJ {
//! particle ids according to PDG
namespace pid {
//! The possible particle identities. We use PDG IDs as standard.
enum ParticleID{
d = 1, /*!< Down Quark */
down = d, /*!< Down Quark */
u = 2, /*!< Up Quark */
up = u, /*!< Up Quark */
s = 3, /*!< Strange Quark */
strange = s, /*!< Strange Quark */
c = 4, /*!< Charm Quark */
charm = c, /*!< Charm Quark */
b = 5, /*!< Bottom Quark */
bottom = b, /*!< Bottom Quark */
t = 6, /*!< Top Quark */
top = t, /*!< Top Quark */
e = 11, /*!< Electron */
electron = e, /*!< Electron */
nu_e = 12, /*!< Electron Neutrino */
electron_neutrino = nu_e, /*!< Electron neutrino */
mu = 13, /*!< Muon */
muon = mu, /*!< Muon */
nu_mu = 14, /*!< Muon Neutrino */
muon_neutrino = nu_mu, /*!< Muon Neutrino */
tau = 15, /*!< Tau */
nu_tau = 16, /*!< Tau Neutrino */
tau_neutrino = nu_tau, /*!< Tau Neutrino */
d_bar = -d, /*!< Anti-Down Quark */
u_bar = -u, /*!< Anti-Up quark */
s_bar = -s, /*!< Anti-Strange Quark */
c_bar = -c, /*!< Anti-Charm Quark */
b_bar = -b, /*!< Anti-Bottom Quark */
t_bar = -t, /*!< Anti-Top Quark */
e_bar = -e, /*!< Positron */
nu_e_bar = -nu_e, /*!< Anti-Electron Neutrino */
mu_bar = -mu, /*!< Anti-Muon */
nu_mu_bar = -nu_mu, /*!< Anti-Muon Neutrino */
tau_bar = -tau, /*!< Anti-Tau */
nu_tau_bar = -nu_tau, /*!< Anti-Tau Neutrino */
gluon = 21, /*!< Gluon */
g = gluon, /*!< Gluon */
photon = 22, /*!< Photon */
gamma = photon, /*!< Photon */
Z = 23, /*!< Z Boson */
Wp = 24, /*!< W- Boson */
Wm = -Wp, /*!< W+ Boson */
h = 25, /*!< Higgs Boson */
Higgs = h, /*!< Higgs Boson */
higgs = h, /*!< Higgs Boson */
p = 2212, /*!< Proton */
proton = p, /*!< Proton */
p_bar = -p, /*!< Anti-Proton */
};
}
using ParticleID = pid::ParticleID;
+ //! Convert a particle name to the corresponding PDG particle ID
+ ParticleID to_ParticleID(std::string const & name);
+
/**
* \brief Function to determine if particle is a parton
* @param p PDG ID of particle
* @returns true if the particle is a parton, false otherwise
*/
inline
constexpr bool is_parton(ParticleID p){
return p == pid::gluon || std::abs(p) <= pid::top;
}
/**
* \brief function to determine if the particle is a photon, W, Z, or Higgs boson
* @param id PDG ID of particle
* @returns true if the partice is a A,W,Z, or H, false otherwise
*/
inline
constexpr bool is_AWZH_boson(ParticleID id){
return id == pid::Wm || (id >= pid::photon && id <= pid::Higgs);
}
}
diff --git a/include/HEJ/YAMLreader.hh b/include/HEJ/YAMLreader.hh
index d2ed33e..8740f14 100644
--- a/include/HEJ/YAMLreader.hh
+++ b/include/HEJ/YAMLreader.hh
@@ -1,244 +1,241 @@
/** \file
* \brief The file which handles the configuration file parameters
*
* The configuration files parameters are read and then stored
* within this objects.
*/
#pragma once
#include "yaml-cpp/yaml.h"
#include "HEJ/config.hh"
#include "HEJ/exceptions.hh"
namespace HEJ{
//! Load configuration from file
/**
* @param config_file Name of the YAML configuration file
* @returns The HEJ 2 configuration
*/
Config load_config(std::string const & config_file);
//! Set option using the corresponding YAML entry
/**
* @param setting Option variable to be set
* @param yaml Root of the YAML configuration
* @param names Name of the entry
*
* If the entry does not exist or has the wrong type or format
* an exception is thrown.
*
* For example
* @code
* set_from_yaml(foobar, yaml, "foo", "bar")
* @endcode
* is equivalent to
* @code
* foobar = yaml["foo"]["bar"].as<decltype(foobar)>()
* @endcode
* with improved diagnostics on errors.
*
* @see set_from_yaml_if_defined
*/
template<typename T, typename... YamlNames>
void set_from_yaml(
T & setting,
YAML::Node const & yaml, YamlNames const & ... names
);
//! Set option using the corresponding YAML entry, if present
/**
* @param setting Option variable to be set
* @param yaml Root of the YAML configuration
* @param names Name of the entry
*
* This function works similar to set_from_yaml, but does not
* throw any exception if the requested YAML entry does not exist.
*
* @see set_from_yaml
*/
template<typename T, typename... YamlNames>
void set_from_yaml_if_defined(
T & setting,
YAML::Node const & yaml, YamlNames const & ... names
);
//! Extract jet parameters from YAML configuration
JetParameters get_jet_parameters(
YAML::Node const & node, std::string const & entry
);
//! Extract Higgs coupling settings from YAML configuration
HiggsCouplingSettings get_Higgs_coupling(
YAML::Node const & node, std::string const & entry
);
//! Extract scale setting parameters from YAML configuration
ScaleConfig to_ScaleConfig(YAML::Node const & yaml);
//! Extract random number generator settings from YAML configuration
RNGConfig to_RNGConfig(YAML::Node const & node, std::string const & entry);
- //! Convert a particle name to the corresponding PDG particle ID
- ParticleID to_ParticleID(std::string const & particle_name);
-
//! Check whether all options in configuration are supported
/**
* @param conf Configuration to be checked
* @param supported Tree of supported options
*
* If conf contains an entry that does not appear in supported
* an unknown_option exception is thrown. Sub-entries of "analysis"
* (if present) are not checked.
*
* @see unknown_option
*/
void assert_all_options_known(
YAML::Node const & conf, YAML::Node const & supported
);
namespace detail{
void set_from_yaml(fastjet::JetAlgorithm & setting, YAML::Node const & yaml);
void set_from_yaml(EventTreatment & setting, YAML::Node const & yaml);
void set_from_yaml(ParticleID & setting, YAML::Node const & yaml);
void set_from_yaml(OutputFile & setting, YAML::Node const & yaml);
inline
void set_from_yaml(YAML::Node & setting, YAML::Node const & yaml){
setting = yaml;
}
template<typename Scalar>
void set_from_yaml(Scalar & setting, YAML::Node const & yaml){
assert(yaml);
if(!yaml.IsScalar()){
throw invalid_type{"value is not a scalar"};
}
try{
setting = yaml.as<Scalar>();
}
catch(...){
throw invalid_type{
"value " + yaml.as<std::string>()
+ " cannot be converted to a " + type_string(setting)
};
}
}
template<typename T>
void set_from_yaml(optional<T> & setting, YAML::Node const & yaml){
T tmp;
set_from_yaml(tmp, yaml);
setting = tmp;
}
template<typename T>
void set_from_yaml(std::vector<T> & setting, YAML::Node const & yaml){
assert(yaml);
// special case: treat a single value like a vector with one element
if(yaml.IsScalar()){
setting.resize(1);
return set_from_yaml(setting.front(), yaml);
}
if(yaml.IsSequence()){
setting.resize(yaml.size());
for(size_t i = 0; i < setting.size(); ++i){
set_from_yaml(setting[i], yaml[i]);
}
return;
}
throw invalid_type{""};
}
template<typename T, typename FirstName, typename... YamlNames>
void set_from_yaml(
T & setting,
YAML::Node const & yaml, FirstName const & name,
YamlNames && ... names
){
if(!yaml[name]) throw missing_option{""};
set_from_yaml(
setting,
yaml[name], std::forward<YamlNames>(names)...
);
}
template<typename T>
void set_from_yaml_if_defined(T & setting, YAML::Node const & yaml){
return set_from_yaml(setting, yaml);
}
template<typename T, typename FirstName, typename... YamlNames>
void set_from_yaml_if_defined(
T & setting,
YAML::Node const & yaml, FirstName const & name,
YamlNames && ... names
){
if(!yaml[name]) return;
set_from_yaml_if_defined(
setting,
yaml[name], std::forward<YamlNames>(names)...
);
}
}
template<typename T, typename... YamlNames>
void set_from_yaml(
T & setting,
YAML::Node const & yaml, YamlNames const & ... names
){
try{
detail::set_from_yaml(setting, yaml, names...);
}
catch(invalid_type const & ex){
throw invalid_type{
"In option " + join(": ", names...) + ": " + ex.what()
};
}
catch(missing_option const &){
throw missing_option{
"No entry for mandatory option " + join(": ", names...)
};
}
catch(std::invalid_argument const & ex){
throw missing_option{
"In option " + join(": ", names...) + ":"
" invalid value " + ex.what()
};
}
}
template<typename T, typename... YamlNames>
void set_from_yaml_if_defined(
T & setting,
YAML::Node const & yaml, YamlNames const & ... names
){
try{
detail::set_from_yaml_if_defined(setting, yaml, names...);
}
catch(invalid_type const & ex){
throw invalid_type{
"In option " + join(": ", names...) + ": " + ex.what()
};
}
catch(std::invalid_argument const & ex){
throw missing_option{
"In option " + join(": ", names...) + ":"
" invalid value " + ex.what()
};
}
}
}
namespace YAML {
template<>
struct convert<HEJ::OutputFile> {
static Node encode(HEJ::OutputFile const & outfile);
static bool decode(Node const & node, HEJ::OutputFile & out);
};
}
diff --git a/src/PDG_codes.cc b/src/PDG_codes.cc
new file mode 100644
index 0000000..3dfc576
--- /dev/null
+++ b/src/PDG_codes.cc
@@ -0,0 +1,48 @@
+#include "HEJ/PDG_codes.hh"
+
+#include <map>
+#include <iostream>
+
+namespace HEJ{
+ ParticleID to_ParticleID(std::string const & name){
+ using namespace HEJ::pid;
+ static const std::map<std::string, ParticleID> known = {
+ {"d", d}, {"down", down}, {"1",static_cast<ParticleID>(1)},
+ {"u", u}, {"up", up}, {"2",static_cast<ParticleID>(2)},
+ {"s", s}, {"strange", strange}, {"3",static_cast<ParticleID>(3)},
+ {"c", c}, {"charm", charm}, {"4",static_cast<ParticleID>(4)},
+ {"b", b}, {"bottom", bottom}, {"5",static_cast<ParticleID>(5)},
+ {"t", t}, {"top", top}, {"6",static_cast<ParticleID>(6)},
+ {"e", e}, {"electron", electron}, {"11",static_cast<ParticleID>(11)},
+ {"nu_e", nu_e}, {"electron_neutrino", electron_neutrino}, {"12",static_cast<ParticleID>(12)},
+ {"mu", mu}, {"muon", muon}, {"13",static_cast<ParticleID>(13)},
+ {"nu_mu", nu_mu}, {"muon_neutrino", muon_neutrino}, {"14",static_cast<ParticleID>(14)},
+ {"tau", tau}, {"15",static_cast<ParticleID>(15)},
+ {"nu_tau", nu_tau}, {"tau_neutrino", tau_neutrino}, {"16",static_cast<ParticleID>(16)},
+ {"d_bar", d_bar}, {"-1",static_cast<ParticleID>(-1)},
+ {"u_bar", u_bar}, {"-2",static_cast<ParticleID>(-2)},
+ {"s_bar", s_bar}, {"-3",static_cast<ParticleID>(-3)},
+ {"c_bar", c_bar}, {"-4",static_cast<ParticleID>(-4)},
+ {"b_bar", b_bar}, {"-5",static_cast<ParticleID>(-5)},
+ {"t_bar", t_bar}, {"-6",static_cast<ParticleID>(-6)},
+ {"e_bar", e_bar}, {"-11",static_cast<ParticleID>(-11)},
+ {"nu_e_bar", nu_e_bar}, {"-12",static_cast<ParticleID>(-12)},
+ {"mu_bar", mu_bar}, {"-13",static_cast<ParticleID>(-13)},
+ {"nu_mu_bar", nu_mu_bar}, {"-14",static_cast<ParticleID>(-14)},
+ {"tau_bar", tau_bar}, {"-15",static_cast<ParticleID>(-15)},
+ {"nu_tau_bar", nu_tau_bar}, {"-16",static_cast<ParticleID>(-16)},
+ {"gluon", gluon}, {"g", g}, {"21",static_cast<ParticleID>(21)},
+ {"photon", photon}, {"gamma", gamma}, {"22",static_cast<ParticleID>(22)},
+ {"Z", Z}, {"23",static_cast<ParticleID>(23)},
+ {"Wp", Wp}, {"W+", Wp}, {"24",static_cast<ParticleID>(24)},
+ {"Wm", Wm}, {"W-", Wm}, {"-24",static_cast<ParticleID>(-24)},
+ {"h", h}, {"H", h}, {"Higgs", Higgs}, {"higgs", higgs}, {"25",static_cast<ParticleID>(25)},
+ {"p", p}, {"proton", proton}, {"p_bar", p_bar}, {"2212",static_cast<ParticleID>(2212)}
+ };
+ const auto res = known.find(name);
+ if(res == known.end()){
+ throw std::invalid_argument("Unknown particle " + name);
+ }
+ return res->second;
+ }
+}
diff --git a/src/YAMLreader.cc b/src/YAMLreader.cc
index ef503de..13f3888 100644
--- a/src/YAMLreader.cc
+++ b/src/YAMLreader.cc
@@ -1,486 +1,462 @@
#include "HEJ/YAMLreader.hh"
#include <set>
#include <string>
#include <vector>
#include <iostream>
#include <stdexcept>
#include <dlfcn.h>
namespace HEJ{
namespace{
//! Get YAML tree of supported options
/**
* The configuration file is checked against this tree of options
* in assert_all_options_known.
*/
YAML::Node const & get_supported_options(){
const static YAML::Node supported = [](){
YAML::Node supported;
static const auto opts = {
"trials", "min extparton pt", "max ext soft pt fraction",
"FKL", "unordered", "non-HEJ",
"scales", "scale factors", "max scale ratio", "import scales",
"log correction", "unweight", "event output", "analysis"
};
// add subnodes to "supported" - the assigned value is irrelevant
for(auto && opt: opts) supported[opt] = "";
for(auto && jet_opt: {"min pt", "algorithm", "R"}){
supported["resummation jets"][jet_opt] = "";
supported["fixed order jets"][jet_opt] = "";
}
for(auto && opt: {"mt", "use impact factors", "include bottom", "mb"}){
supported["Higgs coupling"][opt] = "";
}
for(auto && opt: {"name", "seed"}){
supported["random generator"][opt] = "";
}
return supported;
}();
return supported;
}
fastjet::JetAlgorithm to_JetAlgorithm(std::string const & algo){
using namespace fastjet;
static const std::map<std::string, fastjet::JetAlgorithm> known = {
{"kt", kt_algorithm},
{"cambridge", cambridge_algorithm},
{"antikt", antikt_algorithm},
{"genkt", genkt_algorithm},
{"cambridge for passive", cambridge_for_passive_algorithm},
{"genkt for passive", genkt_for_passive_algorithm},
{"ee kt", ee_kt_algorithm},
{"ee genkt", ee_genkt_algorithm},
{"plugin", plugin_algorithm}
};
const auto res = known.find(algo);
if(res == known.end()){
throw std::invalid_argument("Unknown jet algorithm " + algo);
}
return res->second;
}
EventTreatment to_EventTreatment(std::string const & name){
static const std::map<std::string, EventTreatment> known = {
{"reweight", EventTreatment::reweight},
{"keep", EventTreatment::keep},
{"discard", EventTreatment::discard}
};
const auto res = known.find(name);
if(res == known.end()){
throw std::invalid_argument("Unknown event treatment " + name);
}
return res->second;
}
} // namespace anonymous
- ParticleID to_ParticleID(std::string const & name){
- using namespace HEJ::pid;
- static const std::map<std::string, ParticleID> known = {
- {"d", d}, {"down", down}, {"u", u}, {"up", up}, {"s", s}, {"strange", strange},
- {"c", c}, {"charm", charm}, {"b", b}, {"bottom", bottom}, {"t", t}, {"top", top},
- {"e", e}, {"electron", electron}, {"nu_e", nu_e}, {"electron_neutrino", electron_neutrino},
- {"mu", mu}, {"muon", muon}, {"nu_mu", nu_mu}, {"muon_neutrino", muon_neutrino},
- {"tau", tau}, {"nu_tau", nu_tau}, {"tau_neutrino", tau_neutrino},
- {"d_bar", d_bar}, {"u_bar", u_bar}, {"s_bar", s_bar}, {"c_bar", c_bar},
- {"b_bar", b_bar}, {"t_bar", t_bar}, {"e_bar", e_bar},
- {"nu_e_bar", nu_e_bar}, {"mu_bar", mu_bar}, {"nu_mu_bar", nu_mu_bar},
- {"tau_bar", tau_bar}, {"nu_tau_bar", nu_tau_bar},
- {"gluon", gluon}, {"g", g}, {"photon", photon}, {"gamma", gamma},
- {"Z", Z}, {"Wp", Wp}, {"Wm", Wm}, {"W+", Wp}, {"W-", Wm},
- {"h", h}, {"H", h}, {"Higgs", Higgs}, {"higgs", higgs},
- {"p", p}, {"proton", proton}, {"p_bar", p_bar}
- };
- const auto res = known.find(name);
- if(res == known.end()){
- throw std::invalid_argument("Unknown particle " + name);
- }
- return res->second;
- }
-
namespace detail{
void set_from_yaml(fastjet::JetAlgorithm & setting, YAML::Node const & yaml){
setting = to_JetAlgorithm(yaml.as<std::string>());
}
void set_from_yaml(EventTreatment & setting, YAML::Node const & yaml){
setting = to_EventTreatment(yaml.as<std::string>());
}
void set_from_yaml(ParticleID & setting, YAML::Node const & yaml){
setting = to_ParticleID(yaml.as<std::string>());
}
} // namespace detail
JetParameters get_jet_parameters(
YAML::Node const & node,
std::string const & entry
){
assert(node);
JetParameters result;
fastjet::JetAlgorithm jet_algo = fastjet::antikt_algorithm;
double R;
set_from_yaml_if_defined(jet_algo, node, entry, "algorithm");
set_from_yaml(R, node, entry, "R");
result.def = fastjet::JetDefinition{jet_algo, R};
set_from_yaml(result.min_pt, node, entry, "min pt");
return result;
}
RNGConfig to_RNGConfig(
YAML::Node const & node,
std::string const & entry
){
assert(node);
RNGConfig result;
set_from_yaml(result.name, node, entry, "name");
set_from_yaml_if_defined(result.seed, node, entry, "seed");
return result;
}
HiggsCouplingSettings get_Higgs_coupling(
YAML::Node const & node,
std::string const & entry
){
assert(node);
static constexpr double mt_max = 2e4;
#ifndef HEJ_BUILD_WITH_QCDLOOP
if(node[entry]){
throw std::invalid_argument{
"Higgs coupling settings require building HEJ 2 "
"with QCDloop support"
};
}
#endif
HiggsCouplingSettings settings;
set_from_yaml_if_defined(settings.mt, node, entry, "mt");
set_from_yaml_if_defined(settings.mb, node, entry, "mb");
set_from_yaml_if_defined(settings.include_bottom, node, entry, "include bottom");
set_from_yaml_if_defined(settings.use_impact_factors, node, entry, "use impact factors");
if(settings.use_impact_factors){
if(settings.mt != std::numeric_limits<double>::infinity()){
throw std::invalid_argument{
"Conflicting settings: "
"impact factors may only be used in the infinite top mass limit"
};
}
}
else{
// huge values of the top mass are numerically unstable
settings.mt = std::min(settings.mt, mt_max);
}
return settings;
}
FileFormat to_FileFormat(std::string const & name){
static const std::map<std::string, FileFormat> known = {
{"Les Houches", FileFormat::Les_Houches},
{"HepMC", FileFormat::HepMC}
};
const auto res = known.find(name);
if(res == known.end()){
throw std::invalid_argument("Unknown file format " + name);
}
return res->second;
}
std::string extract_suffix(std::string const & filename){
size_t separator = filename.rfind('.');
if(separator == filename.npos) return {};
return filename.substr(separator + 1);
}
FileFormat format_from_suffix(std::string const & filename){
const std::string suffix = extract_suffix(filename);
if(suffix == "lhe") return FileFormat::Les_Houches;
if(suffix == "hepmc") return FileFormat::HepMC;
throw std::invalid_argument{
"Can't determine format for output file " + filename
};
}
void assert_all_options_known(
YAML::Node const & conf, YAML::Node const & supported
){
if(!conf.IsMap()) return;
if(!supported.IsMap()) throw invalid_type{"must not have sub-entries"};
for(auto const & entry: conf){
const auto name = entry.first.as<std::string>();
if(! supported[name]) throw unknown_option{name};
/* check sub-options, e.g. 'resummation jets: min pt'
* we don't check analysis sub-options
* those depend on the analysis being used and should be checked there
* similar for "import scales"
*/
if(name != "analysis" && name != "import scales"){
try{
assert_all_options_known(conf[name], supported[name]);
}
catch(unknown_option const & ex){
throw unknown_option{name + ": " + ex.what()};
}
catch(invalid_type const & ex){
throw invalid_type{name + ": " + ex.what()};
}
}
}
}
} // namespace HEJ
namespace YAML {
Node convert<HEJ::OutputFile>::encode(HEJ::OutputFile const & outfile) {
Node node;
node[to_string(outfile.format)] = outfile.name;
return node;
};
bool convert<HEJ::OutputFile>::decode(Node const & node, HEJ::OutputFile & out) {
switch(node.Type()){
case NodeType::Map: {
YAML::const_iterator it = node.begin();
out.format = HEJ::to_FileFormat(it->first.as<std::string>());
out.name = it->second.as<std::string>();
return true;
}
case NodeType::Scalar:
out.name = node.as<std::string>();
out.format = HEJ::format_from_suffix(out.name);
return true;
default:
return false;
}
}
} // namespace YAML
namespace HEJ{
namespace detail{
void set_from_yaml(OutputFile & setting, YAML::Node const & yaml){
setting = yaml.as<OutputFile>();
}
}
namespace{
void update_fixed_order_jet_parameters(
JetParameters & fixed_order_jets, YAML::Node const & yaml
){
if(!yaml["fixed order jets"]) return;
set_from_yaml_if_defined(
fixed_order_jets.min_pt, yaml, "fixed order jets", "min pt"
);
fastjet::JetAlgorithm algo = fixed_order_jets.def.jet_algorithm();
set_from_yaml_if_defined(algo, yaml, "fixed order jets", "algorithm");
double R = fixed_order_jets.def.R();
set_from_yaml_if_defined(R, yaml, "fixed order jets", "R");
fixed_order_jets.def = fastjet::JetDefinition{algo, R};
}
// like std::stod, but throw if not the whole string can be converted
double to_double(std::string const & str){
std::size_t pos;
const double result = std::stod(str, &pos);
if(pos < str.size()){
throw std::invalid_argument(str + " is not a valid double value");
}
return result;
}
using EventScale = double (*)(Event const &);
void import_scale_functions(
std::string const & file,
std::vector<std::string> const & scale_names,
std::unordered_map<std::string, EventScale> & known
) {
auto handle = dlopen(file.c_str(), RTLD_NOW);
char * error = dlerror();
if(error != nullptr) throw std::runtime_error{error};
for(auto const & scale: scale_names) {
void * sym = dlsym(handle, scale.c_str());
error = dlerror();
if(error != nullptr) throw std::runtime_error{error};
known.emplace(scale, reinterpret_cast<EventScale>(sym));
}
}
auto get_scale_map(
YAML::Node const & yaml
) {
std::unordered_map<std::string, EventScale> scale_map;
scale_map.emplace("H_T", H_T);
scale_map.emplace("max jet pperp", max_jet_pt);
scale_map.emplace("jet invariant mass", jet_invariant_mass);
scale_map.emplace("m_j1j2", m_j1j2);
if(yaml["import scales"]) {
if(! yaml["import scales"].IsMap()) {
throw invalid_type{"Entry 'import scales' is not a map"};
}
for(auto const & import: yaml["import scales"]) {
const auto file = import.first.as<std::string>();
const auto scale_names =
import.second.IsSequence()
?import.second.as<std::vector<std::string>>()
:std::vector<std::string>{import.second.as<std::string>()};
import_scale_functions(file, scale_names, scale_map);
}
}
return scale_map;
}
// simple (as in non-composite) scale functions
/**
* An example for a simple scale function would be H_T,
* H_T/2 is then composite (take H_T and then divide by 2)
*/
ScaleFunction parse_simple_ScaleFunction(
std::string const & scale_fun,
std::unordered_map<std::string, EventScale> const & known
) {
assert(
scale_fun.empty() ||
(!std::isspace(scale_fun.front()) && !std::isspace(scale_fun.back()))
);
const auto it = known.find(scale_fun);
if(it != end(known)) return {it->first, it->second};
try{
const double scale = to_double(scale_fun);
return {scale_fun, FixedScale{scale}};
} catch(std::invalid_argument const &){}
throw std::invalid_argument{"Unknown scale choice: " + scale_fun};
}
std::string trim_front(std::string const & str){
const auto new_begin = std::find_if(
begin(str), end(str), [](char c){ return ! std::isspace(c); }
);
return std::string(new_begin, end(str));
}
std::string trim_back(std::string str){
size_t pos = str.size() - 1;
// use guaranteed wrap-around behaviour to check whether we have
// traversed the whole string
for(; pos < str.size() && std::isspace(str[pos]); --pos) {}
str.resize(pos + 1); // note that pos + 1 can be 0
return str;
}
ScaleFunction parse_ScaleFunction(
std::string const & scale_fun,
std::unordered_map<std::string, EventScale> const & known
){
assert(
scale_fun.empty() ||
(!std::isspace(scale_fun.front()) && !std::isspace(scale_fun.back()))
);
const size_t delim = scale_fun.find_first_of("*/");
if(delim == scale_fun.npos){
return parse_simple_ScaleFunction(scale_fun, known);
}
const std::string first = trim_back(std::string{scale_fun, 0, delim});
const std::string second = trim_front(std::string{scale_fun, delim+1});
if(scale_fun[delim] == '/'){
return parse_simple_ScaleFunction(first, known)/to_double(second);
}
else{
assert(scale_fun[delim] == '*');
try{
const double factor = to_double(second);
return factor*parse_simple_ScaleFunction(first, known);
}
catch(std::invalid_argument const &){
const double factor = to_double(first);
return factor*parse_simple_ScaleFunction(first, known);
}
}
}
EventTreatMap get_event_treatment(
YAML::Node const & yaml
){
using namespace event_type;
EventTreatMap treat {
{no_2_jets, EventTreatment::discard},
{bad_final_state, EventTreatment::discard},
{FKL, EventTreatment::reweight},
{unob, EventTreatment::keep},
{unof, EventTreatment::keep},
{nonHEJ, EventTreatment::keep}
};
set_from_yaml(treat.at(FKL), yaml, "FKL");
set_from_yaml(treat.at(unob), yaml, "unordered");
treat.at(unof) = treat.at(unob);
set_from_yaml(treat.at(nonHEJ), yaml, "non-HEJ");
if(treat[nonHEJ] == EventTreatment::reweight){
throw std::invalid_argument{"Cannot reweight non-HEJ events"};
}
return treat;
}
Config to_Config(YAML::Node const & yaml){
try{
assert_all_options_known(yaml, get_supported_options());
}
catch(unknown_option const & ex){
throw unknown_option{std::string{"Unknown option '"} + ex.what() + "'"};
}
Config config;
config.resummation_jets = get_jet_parameters(yaml, "resummation jets");
config.fixed_order_jets = config.resummation_jets;
update_fixed_order_jet_parameters(config.fixed_order_jets, yaml);
set_from_yaml(config.min_extparton_pt, yaml, "min extparton pt");
config.max_ext_soft_pt_fraction = std::numeric_limits<double>::infinity();
set_from_yaml_if_defined(
config.max_ext_soft_pt_fraction, yaml, "max ext soft pt fraction"
);
set_from_yaml(config.trials, yaml, "trials");
set_from_yaml(config.log_correction, yaml, "log correction");
set_from_yaml(config.unweight, yaml, "unweight");
config.treat = get_event_treatment(yaml);
set_from_yaml_if_defined(config.output, yaml, "event output");
config.rng = to_RNGConfig(yaml, "random generator");
set_from_yaml_if_defined(config.analysis_parameters, yaml, "analysis");
config.scales = to_ScaleConfig(yaml);
config.Higgs_coupling = get_Higgs_coupling(yaml, "Higgs coupling");
return config;
}
} // namespace anonymous
ScaleConfig to_ScaleConfig(YAML::Node const & yaml){
ScaleConfig config;
auto scale_funs = get_scale_map(yaml);
std::vector<std::string> scales;
set_from_yaml(scales, yaml, "scales");
config.base.reserve(scales.size());
std::transform(
begin(scales), end(scales), std::back_inserter(config.base),
[scale_funs](auto const & entry){
return parse_ScaleFunction(entry, scale_funs);
}
);
set_from_yaml_if_defined(config.factors, yaml, "scale factors");
config.max_ratio = std::numeric_limits<double>::infinity();
set_from_yaml_if_defined(config.max_ratio, yaml, "max scale ratio");
return config;
}
Config load_config(std::string const & config_file){
try{
return to_Config(YAML::LoadFile(config_file));
}
catch(...){
std::cerr << "Error reading " << config_file << ":\n ";
throw;
}
}
} // namespace HEJ

File Metadata

Mime Type
text/x-diff
Expires
Thu, Apr 24, 6:34 AM (1 d, 12 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
4887755
Default Alt Text
(42 KB)

Event Timeline