diff --git a/t/ME_data/config_mt.yml b/t/ME_data/config_mt.yml index 25be9fd..5e1be8d 100644 --- a/t/ME_data/config_mt.yml +++ b/t/ME_data/config_mt.yml @@ -1,44 +1,44 @@ trials: 1 -min extparton pt: 30 +max ext soft pt fraction: 0.1 resummation jets: min pt: 30 algorithm: antikt R: 0.4 fixed order jets: min pt: 30 event treatment: FKL: reweight unordered: reweight extremal qqx: discard central qqx: discard non-resummable: discard scales: 125 log correction: false random generator: name: mixmax seed: 1 vev: 246.2196508 particle properties: Higgs: mass: 125 width: 0.004165 W: mass: 80.385 width: 2.085 Z: mass: 91.187 width: 2.495 Higgs coupling: use impact factors: false mt: 174 include bottom: false diff --git a/t/ME_data/config_mtinf.yml b/t/ME_data/config_mtinf.yml index f1a8f7b..3f87d26 100644 --- a/t/ME_data/config_mtinf.yml +++ b/t/ME_data/config_mtinf.yml @@ -1,39 +1,39 @@ trials: 1 -min extparton pt: 30 +max ext soft pt fraction: 0.1 resummation jets: min pt: 30 algorithm: antikt R: 0.4 fixed order jets: min pt: 30 event treatment: FKL: reweight unordered: reweight extremal qqx: discard central qqx: discard non-resummable: discard scales: 125 log correction: false vev: 246.2196508 particle properties: Higgs: mass: 125 width: 0.004165 W: mass: 80.385 width: 2.085 Z: mass: 91.187 width: 2.495 random generator: name: mixmax seed: 1 diff --git a/t/ME_data/config_mtmb.yml b/t/ME_data/config_mtmb.yml index eb0075a..924f524 100644 --- a/t/ME_data/config_mtmb.yml +++ b/t/ME_data/config_mtmb.yml @@ -1,45 +1,45 @@ trials: 1 -min extparton pt: 30 +max ext soft pt fraction: 0.1 resummation jets: min pt: 30 algorithm: antikt R: 0.4 fixed order jets: min pt: 30 event treatment: FKL: reweight unordered: reweight extremal qqx: discard central qqx: discard non-resummable: discard scales: 125 log correction: false random generator: name: mixmax seed: 1 vev: 246.2196508 particle properties: Higgs: mass: 125 width: 0.004165 W: mass: 80.385 width: 2.085 Z: mass: 91.187 width: 2.495 Higgs coupling: use impact factors: false mt: 174 include bottom: true mb: 4.7 diff --git a/t/ME_data/config_w_ME.yml b/t/ME_data/config_w_ME.yml index 27d826b..d0da43a 100644 --- a/t/ME_data/config_w_ME.yml +++ b/t/ME_data/config_w_ME.yml @@ -1,39 +1,39 @@ trials: 1 -min extparton pt: 30 +max ext soft pt fraction: 0.1 resummation jets: min pt: 30 algorithm: antikt R: 0.4 fixed order jets: min pt: 30 event treatment: FKL: reweight unordered: reweight extremal qqx: reweight central qqx: reweight non-resummable: discard scales: 125 log correction: false vev: 246.2196508 particle properties: Higgs: mass: 125 width: 0.004165 W: mass: 80.419 width: 2.0476 Z: mass: 91.187 width: 2.495 random generator: name: mixmax seed: 1 diff --git a/t/check_res.cc b/t/check_res.cc index d33e50b..e3b8fad 100644 --- a/t/check_res.cc +++ b/t/check_res.cc @@ -1,171 +1,169 @@ /** * \authors The HEJ collaboration (see AUTHORS for details) * \date 2019 * \copyright GPLv2 or later */ #include <cmath> #include <iostream> #include "LHEF/LHEF.h" #include "HEJ/Event.hh" #include "HEJ/EventReweighter.hh" #include "HEJ/Mixmax.hh" #include "HEJ/stream.hh" #define ASSERT(x) if(!(x)) { \ std::cerr << "Assertion '" #x "' failed.\n"; \ return EXIT_FAILURE; \ } namespace{ const fastjet::JetDefinition jet_def{fastjet::kt_algorithm, 0.4}; const fastjet::JetDefinition Born_jet_def{jet_def}; constexpr double Born_jetptmin = 30; - constexpr double extpartonptmin = 30; constexpr double max_ext_soft_pt_fraction = 0.1; constexpr double jetptmin = 35; constexpr bool log_corr = false; const HEJ::ParticleProperties Wprop{80.385, 2.085}; const HEJ::ParticleProperties Zprop{91.187, 2.495}; const HEJ::ParticleProperties Hprop{125, 0.004165}; constexpr double vev = 246.2196508; using EventTreatment = HEJ::EventTreatment; using namespace HEJ::event_type; HEJ::EventTreatMap treat{ {no_2_jets, EventTreatment::discard}, {bad_final_state, EventTreatment::discard}, {non_resummable, EventTreatment::discard}, {unof, EventTreatment::discard}, {unob, EventTreatment::discard}, {qqxexb, EventTreatment::discard}, {qqxexf, EventTreatment::discard}, {qqxmid, EventTreatment::discard}, {FKL, EventTreatment::reweight} }; /// true if colour is allowed for particle bool correct_colour(HEJ::Particle const & part){ if(HEJ::is_AWZH_boson(part) && !part.colour) return true; if(!part.colour) return false; int const colour = part.colour->first; int const anti_colour = part.colour->second; if(part.type == HEJ::ParticleID::gluon) return colour != anti_colour && colour > 0 && anti_colour > 0; if(HEJ::is_quark(part)) return anti_colour == 0 && colour > 0; return colour == 0 && anti_colour > 0; } bool correct_colour(HEJ::Event const & ev){ if(!HEJ::event_type::is_resummable(ev.type())) return true; for(auto const & part: ev.incoming()){ if(!correct_colour(part)) return false; } for(auto const & part: ev.outgoing()){ if(!correct_colour(part)) return false; } return true; } }; int main(int argn, char** argv) { if(argn == 5 && std::string(argv[4]) == "unof"){ --argn; treat[unof] = EventTreatment::reweight; treat[unob] = EventTreatment::discard; treat[FKL] = EventTreatment::discard; } if(argn == 5 && std::string(argv[4]) == "unob"){ --argn; treat[unof] = EventTreatment::discard; treat[unob] = EventTreatment::reweight; treat[FKL] = EventTreatment::discard; } else if(argn == 5 && std::string(argv[4]) == "splitf"){ --argn; treat[qqxexb] = EventTreatment::discard; treat[qqxexf] = EventTreatment::reweight; treat[FKL] = EventTreatment::discard; } else if(argn == 5 && std::string(argv[4]) == "splitb"){ --argn; treat[qqxexb] = EventTreatment::reweight; treat[qqxexf] = EventTreatment::discard; treat[FKL] = EventTreatment::discard; } else if(argn == 5 && std::string(argv[4]) == "qqxmid"){ --argn; treat[qqxmid] = EventTreatment::reweight; treat[FKL] = EventTreatment::discard; } if(argn != 4){ std::cerr << "Usage: check_res eventfile xsection tolerance [uno]"; return EXIT_FAILURE; } const double xsec_ref = std::stod(argv[2]); const double tolerance = std::stod(argv[3]); HEJ::istream in{argv[1]}; LHEF::Reader reader{in}; HEJ::PhaseSpacePointConfig psp_conf; psp_conf.jet_param = HEJ::JetParameters{jet_def, jetptmin}; - psp_conf.min_extparton_pt = extpartonptmin; psp_conf.max_ext_soft_pt_fraction = max_ext_soft_pt_fraction; HEJ::MatrixElementConfig ME_conf; ME_conf.log_correction = log_corr; ME_conf.Higgs_coupling = HEJ::HiggsCouplingSettings{}; ME_conf.ew_parameters.set_vevWZH(vev, Wprop, Zprop, Hprop); HEJ::EventReweighterConfig conf; conf.psp_config = std::move(psp_conf); conf.ME_config = std::move(ME_conf); conf.jet_param = psp_conf.jet_param; conf.treat = treat; reader.readEvent(); const bool has_Higgs = std::find( begin(reader.hepeup.IDUP), end(reader.hepeup.IDUP), 25 ) != end(reader.hepeup.IDUP); const double mu = has_Higgs?125.:91.188; HEJ::ScaleGenerator scale_gen{ {{std::to_string(mu), HEJ::FixedScale{mu}}}, {}, 1. }; HEJ::Mixmax ran{}; HEJ::EventReweighter hej{reader.heprup, std::move(scale_gen), conf, ran}; double xsec = 0.; double xsec_err = 0.; do{ auto ev_data = HEJ::Event::EventData{reader.hepeup}; ev_data.reconstruct_intermediate(); HEJ::Event ev{ ev_data.cluster( Born_jet_def, Born_jetptmin ) }; auto resummed_events = hej.reweight(ev, 100); for(auto const & ev: resummed_events) { ASSERT(correct_colour(ev)); ASSERT(std::isfinite(ev.central().weight)); xsec += ev.central().weight; xsec_err += ev.central().weight*ev.central().weight; } } while(reader.readEvent()); xsec_err = std::sqrt(xsec_err); const double significance = std::abs(xsec - xsec_ref) / std::sqrt( xsec_err*xsec_err + tolerance*tolerance ); std::cout << xsec_ref << " +/- " << tolerance << " ~ " << xsec << " +- " << xsec_err << " => " << significance << " sigma\n"; if(significance > 3.){ std::cerr << "Cross section is off by over 3 sigma!\n"; return EXIT_FAILURE; } return EXIT_SUCCESS; } diff --git a/t/jet_config.yml b/t/jet_config.yml index 9d7dc31..79e0363 100644 --- a/t/jet_config.yml +++ b/t/jet_config.yml @@ -1,41 +1,41 @@ trials: 10 -min extparton pt: 30 +max ext soft pt fraction: 0.1 resummation jets: min pt: 35 algorithm: antikt R: 0.4 fixed order jets: min pt: 30 event treatment: FKL: reweight unordered: keep extremal qqx: discard central qqx: keep non-resummable: keep log correction: false scales: 91.188 random generator: name: ranlux64 vev: 246.2196508 particle properties: Higgs: mass: 125 width: 0.004165 W: mass: 80.385 width: 2.085 Z: mass: 91.187 width: 2.495 event output: - tst.lhe diff --git a/t/jet_config_with_import.yml b/t/jet_config_with_import.yml index 446a58c..027a2e8 100644 --- a/t/jet_config_with_import.yml +++ b/t/jet_config_with_import.yml @@ -1,44 +1,44 @@ trials: 10 -min extparton pt: 30 +max ext soft pt fraction: 0.1 resummation jets: min pt: 35 algorithm: antikt R: 0.4 fixed order jets: min pt: 30 event treatment: FKL: reweight unordered: discard extremal qqx: discard central qqx: discard non-resummable: discard log correction: false scales: softest_jet_pt event output: - tst.lhe random generator: name: ranlux64 vev: 246.2196508 particle properties: Higgs: mass: 125 width: 0.004165 W: mass: 80.385 width: 2.085 Z: mass: 91.187 width: 2.495 import scales: ./libscales.so: softest_jet_pt diff --git a/t/test_psp.cc b/t/test_psp.cc index 92b6a2b..105bc49 100644 --- a/t/test_psp.cc +++ b/t/test_psp.cc @@ -1,68 +1,65 @@ /** * \authors The HEJ collaboration (see AUTHORS for details) * \date 2019 * \copyright GPLv2 or later */ #include "LHEF/LHEF.h" #include "HEJ/stream.hh" #include "HEJ/config.hh" #include "HEJ/event_types.hh" #include "HEJ/Event.hh" #include "HEJ/PhaseSpacePoint.hh" #include "HEJ/Ranlux64.hh" namespace{ constexpr int max_trials = 100; - constexpr double extpartonptmin = 45.; - constexpr double max_ext_soft_pt_fraction = - std::numeric_limits<double>::infinity(); + constexpr double max_ext_soft_pt_fraction = 0.1; const fastjet::JetDefinition jet_def{fastjet::kt_algorithm, 0.4}; constexpr double min_jet_pt = 50; }; int main(int argn, char** argv) { if(argn != 2){ std::cerr << "Usage: " << argv[0] << " eventfile"; return EXIT_FAILURE; } HEJ::istream in{argv[1]}; LHEF::Reader reader{in}; LHEF::Writer writer{std::cerr}; writer.heprup = reader.heprup; HEJ::PhaseSpacePointConfig conf; conf.jet_param = HEJ::JetParameters{jet_def, min_jet_pt}; - conf.min_extparton_pt = extpartonptmin; conf.max_ext_soft_pt_fraction = max_ext_soft_pt_fraction; HEJ::Ranlux64 ran{}; while(reader.readEvent()){ const HEJ::Event ev{ HEJ::Event::EventData{reader.hepeup}( jet_def, min_jet_pt ) }; for(int trial = 0; trial < max_trials; ++trial){ HEJ::PhaseSpacePoint psp{ev, conf, ran}; if(psp.weight() != 0){ HEJ::Event::EventData tmp_ev; tmp_ev.incoming = psp.incoming(); tmp_ev.outgoing = psp.outgoing(); tmp_ev.parameters.central = {0,0,0}; HEJ::Event out_ev{ tmp_ev(jet_def, min_jet_pt) }; if(out_ev.type() != ev.type()){ using HEJ::event_type::name; std::cerr << "Wrong class of phase space point:\n" "original event of class " << name(ev.type()) << ":\n"; writer.hepeup = reader.hepeup; writer.writeEvent(); std::cerr << "Phase space point of class " << name(out_ev.type()) << ":\n"; writer.hepeup = to_HEPEUP(out_ev, &writer.heprup); writer.writeEvent(); return EXIT_FAILURE; } } } } return EXIT_SUCCESS; }