Page MenuHomeHEPForge
Files F9501651

No OneTemporary
Actions

View Options

diff --git a/FixedOrderGen/src/config.cc b/FixedOrderGen/src/config.cc
index 2e3fea7..3cdacfe 100644
--- a/FixedOrderGen/src/config.cc
+++ b/FixedOrderGen/src/config.cc
@@ -1,441 +1,419 @@
/**
* \authors The HEJ collaboration (see AUTHORS for details)
* \date 2019-2020
* \copyright GPLv2 or later
*/
#include "config.hh"
#include <algorithm>
#include <cassert>
#include <cctype>
#include <cmath>
#include <cstdlib>
#include <iostream>
#include <iterator>
#include <stdexcept>
#include "HEJ/PDG_codes.hh"
#include "HEJ/YAMLreader.hh"
#include "HEJ/exceptions.hh"
namespace HEJFOG {
using HEJ::set_from_yaml;
using HEJ::set_from_yaml_if_defined;
using HEJ::pid::ParticleID;
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", "subleading fraction","subleading channels",
"scales", "scale factors", "max scale ratio", "pdf", "vev",
"event output", "analyses", "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", "W", "Z"}){
for(auto && particle_opt: {"mass", "width"}){
supported["particle properties"][particle_type][particle_opt] = "";
}
}
for(auto && particle_type: {"Higgs", "Wp", "W+", "Wm", "W-", "Z"}){
for(auto && particle_opt: {"into", "branching ratio"}){
supported["decays"][particle_type][particle_opt] = "";
}
}
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");
if(result.peak_pt && *result.peak_pt <= result.min_pt)
throw std::invalid_argument{
"Value of option 'peak pt' has to be larger than 'min 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 = 0, end = 0; end != std::string::npos;){
begin = str.find_first_not_of(delims, end);
if(begin == std::string::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', 'H', 'Wm', 'Wp', 'e-', 'e+', 'nu_e', 'nu_e_bar'"
};
}
HEJ::ParticleID reconstruct_boson_id(
std::vector<HEJ::ParticleID> const & ids
){
assert(ids.size()==2);
const int pidsum = ids[0] + ids[1];
if(pidsum == +1) {
assert(HEJ::is_antilepton(ids[0]));
if(HEJ::is_antineutrino(ids[0])) {
throw HEJ::not_implemented{"lepton-flavour violating final state"};
}
assert(HEJ::is_neutrino(ids[1]));
// charged antilepton + neutrino means we had a W+
return HEJ::pid::Wp;
}
if(pidsum == -1) {
assert(HEJ::is_antilepton(ids[0]));
if(HEJ::is_neutrino(ids[1])) {
throw HEJ::not_implemented{"lepton-flavour violating final state"};
}
assert(HEJ::is_antineutrino(ids[0]));
// charged lepton + antineutrino means we had a W-
return HEJ::pid::Wm;
}
throw HEJ::not_implemented{
"final state with leptons "+name(ids[0])+" and "+name(ids[1])
+" not supported"
};
}
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())
&& particles[i].back() == 'j')
result.njets += std::stoi(particles[i]);
else{
const auto pid = HEJ::to_ParticleID(particles[i]);
if(pid==HEJ::pid::Higgs || pid==HEJ::pid::Wp || pid==HEJ::pid::Wm){
if(result.boson)
throw std::invalid_argument{
"More than one outgoing boson is not supported"
};
if(!result.boson_decay.empty())
throw std::invalid_argument{
"Production of a boson together with a lepton is not supported"
};
result.boson = pid;
} else if (HEJ::is_anylepton(pid)){
// Do not accept more leptons, if two leptons are already mentioned
if( result.boson_decay.size()>=2 )
throw std::invalid_argument{"Too many leptons provided"};
if(result.boson)
throw std::invalid_argument{
"Production of a lepton together with a boson is not supported"
};
result.boson_decay.emplace_back(pid);
} else {
throw invalid_outgoing(particles[i]);
}
}
}
if(result.njets < 2){
throw std::invalid_argument{
"Process has to include at least two jets ('j')"
};
}
if(!result.boson_decay.empty()){
std::sort(begin(result.boson_decay),end(result.boson_decay));
assert(!result.boson);
result.boson = reconstruct_boson_id(result.boson_decay);
}
return result;
}
HEJFOG::Subleading to_subleading_channel(YAML::Node const & yaml){
std::string name;
using namespace HEJFOG::subleading;
set_from_yaml(name, yaml);
if(name == "none")
return NONE;
if(name == "all")
return ALL;
Subleading channel = NONE;
if(name == "unordered" || name == "uno")
return channel.set(uno);
if(name == "central qqx" || name == "cqqx")
return channel.set(cqqx);
if(name == "extremal qqx" || name == "eqqx")
return channel.set(eqqx);
throw HEJ::unknown_option("Unknown subleading channel '"+name+"'");
}
Subleading get_subleading_channels(YAML::Node const & node){
using YAML::NodeType;
using namespace HEJFOG::subleading;
// all channels allowed by default
if(!node) return ALL;
switch(node.Type()){
case NodeType::Undefined:
return ALL;
case NodeType::Null:
return NONE;
case NodeType::Scalar:
return to_subleading_channel(node);
case NodeType::Map:
throw HEJ::invalid_type{"map is not a valid option for subleading channels"};
case NodeType::Sequence:
Subleading channels;
for(auto && channel_node: node){
channels |= get_subleading_channels(channel_node);
}
return channels;
}
throw std::logic_error{"unreachable"};
}
Decay get_decay(YAML::Node const & node,
std::string const & entry, std::string const & boson
){
Decay decay;
set_from_yaml(decay.products, node, entry, boson, "into");
decay.branching_ratio=1;
set_from_yaml_if_defined(decay.branching_ratio, node, entry, boson,
"branching ratio");
return decay;
}
std::vector<Decay> get_decays(YAML::Node const & node,
std::string const & entry, std::string const & boson
){
using YAML::NodeType;
if(!node[entry][boson]) return {};
switch(node[entry][boson].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, entry, boson)};
case NodeType::Sequence:
std::vector<Decay> result;
for(auto && decay_str: node[entry][boson]){
result.emplace_back(get_decay(decay_str, entry, boson));
}
return result;
}
throw std::logic_error{"unreachable"};
}
ParticlesDecayMap get_all_decays(YAML::Node const & node,
std::string const & entry
){
if(!node[entry]) return {};
if(!node[entry].IsMap())
throw HEJ::invalid_type{entry + " have to be a map"};
ParticlesDecayMap result;
for(auto const & sub_node: node[entry]) {
const auto boson = sub_node.first.as<std::string>();
const auto id = HEJ::to_ParticleID(boson);
result.emplace(id, get_decays(node, entry, boson));
}
return result;
}
- HEJ::ParticleProperties get_particle_properties(
- YAML::Node const & node, std::string const & entry,
- std::string const & boson
- ){
- HEJ::ParticleProperties result{};
- set_from_yaml(result.mass, node, entry, boson, "mass");
- set_from_yaml(result.width, node, entry, boson, "width");
- return result;
- }
-
- HEJ::EWConstants get_ew_parameters(YAML::Node const & node){
- HEJ::EWConstants result;
- double vev = NAN;
- set_from_yaml(vev, node, "vev");
- result.set_vevWZH(vev,
- get_particle_properties(node, "particle properties", "W"),
- get_particle_properties(node, "particle properties", "Z"),
- get_particle_properties(node, "particle properties", "Higgs")
- );
- 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){
auto const & 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.subleading_fraction, yaml, "subleading fraction");
if(config.subleading_fraction == 0)
config.subleading_channels.reset();
else
config.subleading_channels = get_subleading_channels(yaml["subleading channels"]);
- config.ew_parameters = get_ew_parameters(yaml);
+ config.ew_parameters = HEJ::get_ew_parameters(yaml);
config.particle_decays = get_all_decays(yaml, "decays");
if(config.process.boson // check that Ws always decay
&& std::abs(*config.process.boson) == HEJ::ParticleID::Wp
&& config.process.boson_decay.empty()
){
auto const & decay = config.particle_decays.find(*config.process.boson);
if(decay == config.particle_decays.end() || decay->second.empty())
throw std::invalid_argument{
"Decay for "+name(*config.process.boson)+" is required"};
}
set_from_yaml_if_defined(config.analyses_parameters, yaml, "analyses");
if(yaml["analysis"].IsDefined()){
std::cerr <<
"WARNING: Configuration entry 'analysis' is deprecated. "
" Use 'analyses' instead.\n";
YAML::Node ana;
set_from_yaml(ana, yaml, "analysis");
if(!ana.IsNull()){
config.analyses_parameters.push_back(ana);
}
}
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"].IsDefined()) config.unweight = get_unweight(yaml, "unweight");
return config;
}
} // namespace
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 HEJFOG

File Metadata

Mime Type
text/x-diff
Expires
Sun, Feb 23, 2:53 PM (1 d, 7 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
4486728
Default Alt Text
(16 KB)

Event Timeline