Page MenuHomeHEPForge

No OneTemporary

This file is larger than 256 KB, so syntax highlighting was skipped.
diff --git a/FixedOrderGen/CMakeLists.txt b/FixedOrderGen/CMakeLists.txt
new file mode 100644
index 0000000..9745702
--- /dev/null
+++ b/FixedOrderGen/CMakeLists.txt
@@ -0,0 +1,82 @@
+cmake_minimum_required(VERSION 2.8 FATAL_ERROR)
+
+set(CMAKE_LEGACY_CYGWIN_WIN32 0)
+set(CMAKE_EXPORT_COMPILE_COMMANDS ON)
+
+project("Fixed Order Generation" C CXX)
+
+## Flags for the compiler. No warning allowed.
+if ("${CMAKE_CXX_COMPILER_ID}" STREQUAL "GNU")
+ set(warnings "-Wall -Wextra -Werror -std=c++1y")
+elseif ("${CMAKE_CXX_COMPILER_ID}" STREQUAL "Clang")
+ set(warnings "-Wall -Wextra -Werror -std=c++11")
+elseif ("${CMAKE_CXX_COMPILER_ID}" STREQUAL "MSVC")
+ set(warnings "/W4 /WX /EHsc")
+endif()
+if (NOT CONFIGURED_ONCE)
+ set(CMAKE_CXX_FLAGS "${warnings}"
+ CACHE STRING "Flags used by the compiler during all build types." FORCE)
+ set(CMAKE_C_FLAGS "${warnings}"
+ CACHE STRING "Flags used by the compiler during all build types." FORCE)
+endif()
+
+
+## Add directories and find dependences
+include_directories(${CMAKE_CURRENT_SOURCE_DIR}/include ${CMAKE_CURRENT_SOURCE_DIR}/../include)
+set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} "${CMAKE_SOURCE_DIR}/../cmake/Modules/")
+
+#find_library(rhej rhej_)
+find_package(ROOT MODULE COMPONENTS Hist Tree MathCore)
+find_package(fastjet REQUIRED)
+find_package(clhep REQUIRED)
+find_package(lhapdf REQUIRED)
+find_package(gsl REQUIRED)
+find_package(Boost REQUIRED COMPONENTS iostreams)
+find_package(HepMC 3)
+if(${HepMC_FOUND})
+ set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -DRHEJ_BUILD_WITH_HepMC")
+endif()
+find_package(QCDloop 2)
+if(${QCDloop_FOUND})
+ set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -DRHEJ_BUILD_WITH_QCDLOOP")
+ include_directories(SYSTEM ${QCDloop_INCLUDE_PATH})
+endif()
+
+include_directories(SYSTEM ${ROOT_INCLUDE_DIR})
+include_directories(SYSTEM ${lhapdf_INCLUDE_PATH})
+include_directories(SYSTEM ${fastjet_INCLUDE_PATH})
+include_directories(SYSTEM ${clhep_INCLUDE_PATH})
+include_directories(SYSTEM ${gsl_INCLUDE_PATH})
+link_directories(${ROOT_LIBRARY_DIR})
+
+#add_subdirectory(src)
+
+## define executable
+file(GLOB FOgen_source ${CMAKE_CURRENT_SOURCE_DIR}/src/*.cc)
+list(REMOVE_ITEM FOgen_source ${CMAKE_CURRENT_SOURCE_DIR}/src/main.cc)
+add_library(hejfog STATIC ${FOgen_source})
+add_executable(FOgen ${CMAKE_CURRENT_SOURCE_DIR}/src/main.cc)
+
+## link libraries
+set(libraries ${CMAKE_DL_LIBS} ${LHAPDF_LIBRARIES} ${CLHEP_LIBRARIES} ${FASTJET_LIBRARIES} ${GSL_LIBRARIES} ${Boost_LIBRARIES} ${HepMC_LIBRARIES} ${ROOT_LIBRARIES} yaml-cpp)
+if(${QCDloop_FOUND})
+ list(APPEND libraries ${QCDloop_LIBRARIES} quadmath)
+endif()
+# add libraries for reversed HEJ <by hand>
+list(APPEND libraries "/home/andersen/HEJ/PURE/GitReverse/run/lib/librhej.so")
+target_link_libraries(hejfog ${libraries})
+target_link_libraries(FOgen hejfog)
+
+install(TARGETS FOgen DESTINATION bin)
+
+enable_testing()
+set(tst_dir "${CMAKE_CURRENT_SOURCE_DIR}/t")
+add_executable(test_h_2j ${tst_dir}/h_2j.cc)
+target_link_libraries(test_h_2j hejfog)
+add_test(NAME t_h_2j COMMAND test_h_2j WORKING_DIRECTORY ${tst_dir})
+add_executable(test_h_3j ${tst_dir}/h_3j.cc)
+target_link_libraries(test_h_3j hejfog)
+add_test(NAME t_h_3j COMMAND test_h_3j WORKING_DIRECTORY ${tst_dir})
+
+set(CONFIGURED_ONCE TRUE CACHE INTERNAL
+ "A flag showing that CMake has configured at least once.")
diff --git a/FixedOrderGen/configFO.yml b/FixedOrderGen/configFO.yml
new file mode 100644
index 0000000..fa3195c
--- /dev/null
+++ b/FixedOrderGen/configFO.yml
@@ -0,0 +1,42 @@
+trials : 200
+
+jets:
+ min pt: 30
+ R: 0.4
+ algorithm: antikt
+ max rapidity: 5
+
+beam:
+ energy: 6500
+ particles: [p, p]
+
+pdf: 230000
+
+process: p p => h 4j
+
+unordered: true
+
+scales: max jet pperp
+
+event output:
+ - HEJ.lhe
+# - RHEJ.lhe
+# - RHEJ_events.hepmc3
+
+# analysis:
+# # to use a custom analysis
+# plugin: /home/andersen/HEJ/PURE/GitReverse/build/analysis-plugins/src/libVBF.so
+# min dy12: 2.8
+# min m12: 400
+# output: HEJFO.root
+# # wtwt cut: # optional cut on (event weight)^2
+
+#RanLux init: ranlux.0 # file for initialisation of random number engine
+
+# parameters for Higgs-gluon couplings
+# this requires compilation with looptools
+# Higgs coupling:
+# use impact factors: false
+# mt: 174
+# include bottom: true
+# mb: 4.7
diff --git a/FixedOrderGen/include/PhaseSpacePoint.hh b/FixedOrderGen/include/PhaseSpacePoint.hh
new file mode 100644
index 0000000..d2b917f
--- /dev/null
+++ b/FixedOrderGen/include/PhaseSpacePoint.hh
@@ -0,0 +1,119 @@
+/** \file PhaseSpacePoint.hh
+ * \brief Contains the PhaseSpacePoint Class
+ */
+
+#pragma once
+
+#include <vector>
+
+#include <CLHEP/Random/Randomize.h>
+#include <CLHEP/Random/RanluxEngine.h>
+
+#include "fastjet/JetDefinition.hh"
+#include "RHEJ/utility.hh"
+#include "RHEJ/Event.hh"
+#include "RHEJ/PDG_codes.hh"
+#include "RHEJ/PDF.hh"
+
+namespace HEJFOG{
+ class Process;
+
+ using RHEJ::Sparticle;
+ //! A point in resummation phase space
+ class PhaseSpacePoint{
+ public:
+ //! Default PhaseSpacePoint Constructor
+ PhaseSpacePoint() = default;
+
+ //! PhaseSpacePoint Constructor
+ /**
+ * @param proc The process to generate
+ * @param jet_def The Jet Definition Used
+ * @param jetptmin Minimum Jet Transverse Momentum
+ * @param rapmax Maximum parton rapidity
+ * @param pdf The pdf set (used for sampling)
+ */
+ PhaseSpacePoint(
+ Process const & proc,
+ fastjet::JetDefinition jet_def,double jetptmin, double rapmax,
+ RHEJ::PDF & pdf, double E_beam
+ );
+
+ //! Get Weight Function
+ /**
+ * @returns Weight of Event
+ */
+ double weight() const{
+ return weight_;
+ }
+
+ double status() const{
+ return status_;
+ }
+
+ //! Get Incoming Function
+ /**
+ * @returns Incoming Particles
+ */
+ std::array<Sparticle, 2> const & incoming() const{
+ return incoming_;
+ }
+
+ //! Get Outgoing Function
+ /**
+ * @returns Outgoing Particles
+ */
+ std::vector<Sparticle> const & outgoing() const{
+ return outgoing_;
+ }
+
+ static void reset_ranlux(std::string const & init_file);
+ static void reset_ranlux(char const * init_file);
+
+ private:
+
+ std::vector<fastjet::PseudoJet> gen_LO_partons(
+ int count, double ptmin, double ptmax, double rapmax
+ );
+ Sparticle gen_boson(
+ RHEJ::ParticleID bosonid
+ );
+ bool jets_ok(
+ std::vector<fastjet::PseudoJet> const & Born_jets,
+ std::vector<fastjet::PseudoJet> const & partons
+ ) const;
+ void reconstruct_incoming(
+ std::array<RHEJ::ParticleID, 2> const & ids,
+ RHEJ::PDF & pdf, double E_beam
+ );
+ RHEJ::ParticleID generate_incoming_id(double x, double uf, RHEJ::PDF & pdf);
+
+ bool momentum_conserved(double ep) const;
+
+ RHEJ::Sparticle const & most_backward_FKL(
+ std::vector<RHEJ::Sparticle> const & partons
+ ) const;
+ RHEJ::Sparticle const & most_forward_FKL(
+ std::vector<RHEJ::Sparticle> const & partons
+ ) const;
+ RHEJ::Sparticle & most_backward_FKL(std::vector<RHEJ::Sparticle> & partons) const;
+ RHEJ::Sparticle & most_forward_FKL(std::vector<RHEJ::Sparticle> & partons) const;
+ bool extremal_FKL_ok(
+ std::vector<fastjet::PseudoJet> const & partons
+ ) const;
+
+ bool unob_, unof_;
+ double weight_;
+
+ int status_;
+ double jetptmin_;
+
+ fastjet::JetDefinition jet_def_;
+
+ std::array<Sparticle, 2> incoming_;
+ std::vector<Sparticle> outgoing_;
+
+ static CLHEP::Ranlux64Engine ran_;
+ };
+ RHEJ::UnclusteredEvent to_UnclusteredEvent(PhaseSpacePoint const & psp);
+}
diff --git a/FixedOrderGen/include/Process.hh b/FixedOrderGen/include/Process.hh
new file mode 100644
index 0000000..29096e8
--- /dev/null
+++ b/FixedOrderGen/include/Process.hh
@@ -0,0 +1,15 @@
+#pragma once
+
+#include <array>
+
+#include "RHEJ/PDG_codes.hh"
+#include "RHEJ/optional.hh"
+
+namespace HEJFOG{
+ struct Process{
+ std::array<RHEJ::ParticleID, 2> incoming;
+ int njets;
+ RHEJ::optional<RHEJ::ParticleID> boson;
+ };
+
+}
diff --git a/FixedOrderGen/include/config.hh b/FixedOrderGen/include/config.hh
new file mode 100644
index 0000000..c52947a
--- /dev/null
+++ b/FixedOrderGen/include/config.hh
@@ -0,0 +1,47 @@
+#pragma once
+
+#include <array>
+
+#include "yaml-cpp/yaml.h"
+
+#include "fastjet/JetDefinition.hh"
+#include "RHEJ/HiggsCouplingSettings.hh"
+#include "RHEJ/PDG_codes.hh"
+#include "RHEJ/optional.hh"
+#include "RHEJ/config.hh"
+#include "RHEJ/output_formats.hh"
+#include "RHEJ/exceptions.hh"
+
+#include "Process.hh"
+
+namespace HEJFOG{
+ struct JetParameters{
+ fastjet::JetDefinition def;
+ double min_pt;
+ double max_y;
+ };
+
+ struct Beam{
+ double energy;
+ std::array<RHEJ::ParticleID, 2> particles = {
+ RHEJ::pid::proton, RHEJ::pid::proton
+ };
+ };
+
+ struct Config{
+ Process process;
+ int trials;
+ JetParameters jets;
+ Beam beam;
+ int pdf_id;
+ bool unordered = true;
+ YAML::Node analysis_parameters;
+ RHEJ::ScaleGenerator scale_gen;
+ std::vector<RHEJ::OutputFile> output;
+ RHEJ::optional<std::string> RanLux_init;
+ RHEJ::HiggsCouplingSettings Higgs_coupling;
+ };
+
+ Config load_config(std::string const & config_file);
+
+}
diff --git a/FixedOrderGen/src/PhaseSpacePoint.cc b/FixedOrderGen/src/PhaseSpacePoint.cc
new file mode 100644
index 0000000..2798edf
--- /dev/null
+++ b/FixedOrderGen/src/PhaseSpacePoint.cc
@@ -0,0 +1,343 @@
+#include "PhaseSpacePoint.hh"
+
+#include <random>
+
+#include "RHEJ/kinematics.hh"
+#include "RHEJ/utility.hh"
+#include "RHEJ/debug.hh"
+
+#include "Process.hh"
+
+using namespace RHEJ;
+
+namespace HEJFOG{
+
+ static_assert(
+ std::numeric_limits<double>::has_quiet_NaN,
+ "no quiet NaN for double"
+ );
+ constexpr double NaN = std::numeric_limits<double>::quiet_NaN();
+ RHEJ::UnclusteredEvent to_UnclusteredEvent(PhaseSpacePoint const & psp){
+ RHEJ::UnclusteredEvent result;
+ result.incoming = psp.incoming();
+ std::sort(
+ begin(result.incoming), end(result.incoming),
+ [](Sparticle o1, Sparticle o2){return o1.p.pz()<o2.p.pz();}
+ );
+ assert(result.incoming.size() == 2);
+ result.outgoing = psp.outgoing();
+ assert(
+ std::is_sorted(
+ begin(result.outgoing), end(result.outgoing),
+ RHEJ::rapidity_less{}
+ )
+ );
+ assert(result.outgoing.size() >= 2);
+ result.central.mur = NaN;
+ result.central.muf = NaN;
+ result.central.weight = psp.weight();
+ return result;
+ }
+
+
+
+namespace{
+ //generate Ranlux64Engine with fixed, predefined state
+ /*
+ * some (all?) of the Ranlux64Engine constructors leave fields
+ * uninitialised, invoking undefined behaviour. This can be
+ * circumvented by restoring the state from a file
+ */
+ CLHEP::Ranlux64Engine gen_Ranlux64Engine(){
+ static const std::string state =
+ "9876\n"
+ "0.91280703978419097666\n"
+ "0.41606065829518357191\n"
+ "0.99156342622341142601\n"
+ "0.030922955274050423213\n"
+ "0.16206278421638486975\n"
+ "0.76151768001958330956\n"
+ "0.43765760066092695979\n"
+ "0.42904698253748563275\n"
+ "0.11476317525663759511\n"
+ "0.026620053590963976831\n"
+ "0.65953715764414511114\n"
+ "0.30136722624439826745\n"
+ "3.5527136788005009294e-15 4\n"
+ "1 202\n";
+ const std::string file = std::tmpnam(nullptr);
+ {
+ std::ofstream out{file};
+ out << state;
+ }
+ CLHEP::Ranlux64Engine result;
+ result.restoreStatus(file.c_str());
+ return result;
+ }
+ }
+
+ CLHEP::Ranlux64Engine PhaseSpacePoint::ran_{gen_Ranlux64Engine()};
+
+
+ void PhaseSpacePoint::reset_ranlux(std::string const & init_file){
+ reset_ranlux(init_file.c_str());
+ }
+
+ void PhaseSpacePoint::reset_ranlux(char const * init_file){
+ ran_.restoreStatus(init_file);
+ }
+
+ PhaseSpacePoint::PhaseSpacePoint(
+ Process const & proc,
+ fastjet::JetDefinition jet_def,double jetptmin, double rapmax,
+ RHEJ::PDF & pdf, double E_beam
+ ):
+ unob_(false),
+ unof_(false),
+ jetptmin_{jetptmin},
+ jet_def_{jet_def}
+ {
+ const int nout = proc.njets + (proc.boson?1:0);
+ status_ = 0;
+ weight_ = 1;
+ weight_ /= std::tgamma(nout);
+ {
+
+ outgoing_.reserve(nout);
+
+ // sqrt(s)/2 is the maximum pt
+ for(auto&& p_out: gen_LO_partons(proc.njets, jetptmin_, E_beam, rapmax)){
+ outgoing_.emplace_back(Sparticle{pid::gluon, std::move(p_out)});
+ }
+ if(status_ != 0) return;
+
+ if(proc.boson && *proc.boson == pid::Higgs){
+ // The Higgs
+ auto Hparticle=gen_boson(pid::higgs);
+ auto pos=std::upper_bound(
+ begin(outgoing_),end(outgoing_),Hparticle,rapidity_less{}
+ );
+ outgoing_.insert(pos,Hparticle);
+ }
+
+ reconstruct_incoming(proc.incoming, pdf, E_beam);
+ if(status_ != 0) return;
+ // set outgoing states
+ most_backward_FKL(outgoing_).type = incoming_[0].type;
+ most_forward_FKL(outgoing_).type = incoming_[1].type;
+ }
+ }
+
+ std::vector<fastjet::PseudoJet> PhaseSpacePoint::gen_LO_partons(
+ int np , double ptmin, double ptmax, double rapmax
+ ){
+ if (np<2) throw std::invalid_argument{"Not enough partons in gen_LO_partons"};
+ // heuristic parameters for pt sampling
+ const double ptpar = ptmin + np/5.;
+ const double temp1 = atan((ptmax - ptmin)/ptpar);
+
+ weight_ /= pow(16.*pow(M_PI,3),np-2);
+ std::vector<fastjet::PseudoJet> partons(np);
+ for(size_t i = 0; i < (size_t) np; ++i){
+ const double r1 = ran_.flat();
+ const double pt = ptmin + ptpar*tan(r1*temp1);
+ const double temp2 = cos(r1*temp1);
+ const double phi = 2*M_PI*ran_.flat();
+ weight_ *= 2.0*M_PI*pt*ptpar*temp1/(temp2*temp2);
+ const double y = -rapmax + 2*rapmax*ran_.flat();
+ weight_ *= 2*rapmax;
+ partons[i].reset_PtYPhiM(pt, y, phi);
+ partons[i].set_user_index(i + 1);
+
+ assert(ptmin-1e-5 <= partons[i].pt() && partons[i].pt() <= ptmax+1e-5);
+ }
+ // Need to check that at LO, the number of jets = number of partons;
+ fastjet::ClusterSequence cs(partons, jet_def_);
+ auto cluster_jets=cs.inclusive_jets(jetptmin_);
+ if (cluster_jets.size()!=unsigned(np)){
+ weight_=0.0;
+ status_ = 2;
+ return {};
+ }
+
+ std::sort(begin(partons), end(partons), rapidity_less{});
+
+ return partons;
+ }
+
+ Sparticle PhaseSpacePoint::gen_boson(
+ RHEJ::ParticleID bosonid
+ ){
+ std::array<double,2> ptrans{0.,0.};
+ for (auto const & parton:outgoing_) {
+ ptrans[0]-=parton.px();
+ ptrans[1]-=parton.py();
+ }
+
+ // The Higgs:
+ // Generate a y Gaussian distributed around 0
+ const double r1=ran_.flat();
+ const double r2=ran_.flat();
+ const double lninvr1=-log(r1);
+ constexpr double a=1.6;
+ const double temp=a*sqrt(2.*lninvr1)*cos(2.*M_PI*r2);
+ const double y=temp;
+ weight_=weight_*(exp(temp*temp/2./a/a))*sqrt(2.*M_PI)*a;
+
+ // r1=ran.flat();
+ // double sH=flags.mH*(flags.mH + flags.GammaH*tan((M_PI*r1)/2. + (-1. + r1)*atan(flags.mH/flags.GammaH)));
+ const double sH=125.*125.;
+
+ const double mHperp=sqrt(ptrans[0]*ptrans[0]+ptrans[1]*ptrans[1]+sH);
+ const double pz=mHperp*sinh(y);
+ const double E=mHperp*cosh(y);
+
+ return Sparticle{bosonid,fastjet::PseudoJet{ptrans[0],ptrans[1],pz,E}};
+ }
+
+ Sparticle const & PhaseSpacePoint::most_backward_FKL(
+ std::vector<Sparticle> const & partons
+ ) const{
+ if(unob_) return partons[1];
+ if(!RHEJ::is_parton(partons[0])) return partons[1];
+ return partons[0];
+ }
+
+ Sparticle const & PhaseSpacePoint::most_forward_FKL(
+ std::vector<Sparticle> const & partons
+ ) const{
+ const size_t last_idx = partons.size() - 1;
+ if(unof_) return partons[last_idx-1];
+ if(!RHEJ::is_parton(partons[last_idx])) return partons[last_idx-1];
+ return partons[last_idx];
+ }
+
+ Sparticle & PhaseSpacePoint::most_backward_FKL(
+ std::vector<Sparticle> & partons
+ ) const{
+ if(unob_) return partons[1];
+ if(!RHEJ::is_parton(partons[0])) return partons[1];
+ return partons[0];
+ }
+
+ Sparticle & PhaseSpacePoint::most_forward_FKL(
+ std::vector<Sparticle> & partons
+ ) const{
+ const size_t last_idx = partons.size() - 1;
+ if(unof_) return partons[last_idx-1];
+ if(!RHEJ::is_parton(partons[last_idx])) return partons[last_idx-1];
+ return partons[last_idx];
+ }
+
+ void PhaseSpacePoint::reconstruct_incoming(
+ std::array<RHEJ::ParticleID, 2> const & ids,
+ RHEJ::PDF & pdf, double E_beam
+ ){
+ std::tie(incoming_[0].p, incoming_[1].p) = incoming_momenta(outgoing_);
+ // calculate xa, xb
+ const double sqrts=2*E_beam;
+ const double xa=(incoming_[0].p.e()-incoming_[0].p.pz())/sqrts;
+ const double xb=(incoming_[1].p.e()+incoming_[1].p.pz())/sqrts;
+ // abort if phase space point is outside of collider energy reach
+ if (xa>1. || xb>1.){
+ weight_=0;
+ status_ = 3;
+ return;
+ }
+ // pick pdfs
+ /** TODO:
+ * ufa, ufb don't correspond to our final scale choice.
+ * The reversed HEJ scale generators currently expect a full event as input,
+ * so fixing this is not completely trivial
+ */
+ if(ids[0] == pid::proton || ids[0] == pid::p_bar){
+ const double ufa=jetptmin_;
+ incoming_[0].type = generate_incoming_id(xa, ufa, pdf);
+ }
+ else {
+ incoming_[0].type = ids[0];
+ }
+ if(ids[1] == pid::proton || ids[1] == pid::p_bar){
+ const double ufb=jetptmin_;
+ incoming_[1].type = generate_incoming_id(xb, ufb, pdf);
+ }
+ else {
+ incoming_[1].type = ids[1];
+ }
+ assert(momentum_conserved(1e-7));
+ }
+
+ RHEJ::ParticleID PhaseSpacePoint::generate_incoming_id(
+ double x, double uf, RHEJ::PDF & pdf
+ ){
+ const double pdfg=fabs(pdf.pdfpt(0,x,uf,pid::gluon));
+ const double pdfu=fabs(pdf.pdfpt(0,x,uf,pid::up));
+ const double pdfd=fabs(pdf.pdfpt(0,x,uf,pid::down));
+ const double pdfux=fabs(pdf.pdfpt(0,x,uf,pid::u_bar));
+ const double pdfdx=fabs(pdf.pdfpt(0,x,uf,pid::d_bar));
+ const double pdfc=fabs(pdf.pdfpt(0,x,uf,pid::charm));
+ const double pdfs=fabs(pdf.pdfpt(0,x,uf,pid::strange));
+ const double pdfsx=fabs(pdf.pdfpt(0,x,uf,pid::s_bar));
+ const double pdfb=fabs(pdf.pdfpt(0,x,uf,pid::b));
+
+ const double pdftot=pdfg+4.0/9.0*(pdfu + pdfd + pdfux + pdfdx +pdfs +pdfsx + 2.0*(pdfc +pdfb ));
+ const double r1=pdftot*ran_.flat();
+ double sum;
+
+ if (r1<(sum=pdfg)) {
+ weight_*=pdftot/pdfg;
+ return pid::gluon;
+ }
+ if (r1<(sum+=(4./9.)*pdfu)) {
+ weight_*=pdftot/((4./9.)*pdfu);
+ return pid::up;
+ }
+ else if (r1<(sum+=(4./9.)*pdfd)) {
+ weight_*=pdftot/((4./9.)*pdfd);
+ return pid::down;
+ }
+ else if (r1<(sum+=(4./9.)*pdfux)) {
+ weight_*=pdftot/((4./9.)*pdfux);
+ return pid::u_bar;
+ }
+ else if (r1<(sum+=(4./9.)*pdfdx)) {
+ weight_*=pdftot/((4./9.)*pdfdx);
+ return pid::d_bar;
+ }
+ else if (r1<(sum+=(4./9.)*pdfc)) {
+ weight_*=pdftot/((4./9.)*pdfc);
+ return pid::c;
+ }
+ else if (r1<(sum+=(4./9.)*pdfc)){
+ weight_*=pdftot/((4./9.)*pdfc);
+ return pid::c_bar;
+ }
+ else if (r1<(sum+=(4./9.)*pdfs)) {
+ weight_*=pdftot/((4./9.)*pdfs);
+ return pid::s;
+ }
+ else if (r1<(sum+=(4./9.)*pdfsx)) {
+ weight_*=pdftot/((4./9.)*pdfsx);
+ return pid::s_bar;
+ }
+ else if (r1<(sum+=(4./9.)*pdfb)) {
+ weight_*=pdftot/((4./9.)*pdfb);
+ return pid::b;
+ }
+ else if (r1<=(sum+=(4./9.)*pdfb)) {
+ weight_*=pdftot/((4./9.)*pdfb);
+ return pid::b_bar;
+ }
+ std::cout << "Error in choosing incoming parton: "<<x<<" "<<uf<<" "<<sum<<" "<<pdftot<<" "<<r1;
+ std::cout << " "<<pdfg+4./9.*(pdfu+pdfux+pdfd+pdfdx+pdfs+pdfsx+2.*(pdfc+pdfb))<<std::endl;
+ throw std::logic_error{"Failed to choose parton flavour"};
+ }
+
+ bool PhaseSpacePoint::momentum_conserved(double ep) const{
+ fastjet::PseudoJet diff;
+ for(auto const & in: incoming()) diff += in.p;
+ for(auto const & out: outgoing()) diff -= out.p;
+ return nearby_ep(diff, fastjet::PseudoJet{}, ep);
+ }
+
+}
diff --git a/FixedOrderGen/src/config.cc b/FixedOrderGen/src/config.cc
new file mode 100644
index 0000000..ef6d1a3
--- /dev/null
+++ b/FixedOrderGen/src/config.cc
@@ -0,0 +1,205 @@
+#include "config.hh"
+
+#include <cctype>
+
+#include "RHEJ/config.hh"
+
+namespace HEJFOG{
+ using RHEJ::set_from_yaml;
+ using RHEJ::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", "trials", "unordered", "scales", "scale factors",
+ "max scale ratio", "pdf", "event output", "analysis", "RanLux init"
+ };
+ // add subnodes to "supported" - the assigned value is irrelevant
+ for(auto && opt: opts) supported[opt] = "";
+ for(auto && jet_opt: {"min pt", "algorithm", "R", "max rapidity"}){
+ supported["jets"][jet_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] = "";
+ }
+ return supported;
+ }();
+ return supported;
+ }
+
+ JetParameters get_jet_parameters(
+ YAML::Node const & node, std::string const & entry
+ ){
+ const auto p = RHEJ::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");
+ return result;
+ }
+
+ Beam get_Beam(
+ YAML::Node const & node, std::string const & entry
+ ){
+ Beam beam;
+ std::vector<RHEJ::ParticleID> particles;
+ set_from_yaml(beam.energy, node, entry, "energy");
+ set_from_yaml_if_defined(particles, node, entry, "particles");
+ if(! particles.empty()){
+ for(RHEJ::ParticleID particle: particles){
+ if(particle != RHEJ::pid::p && particle != RHEJ::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() = RHEJ::to_ParticleID(particles[0]);
+ result.incoming.back() = RHEJ::to_ParticleID(particles[1]);
+
+ for(size_t i = 0; i < result.incoming.size(); ++i){
+ const RHEJ::ParticleID in = result.incoming[i];
+ if(
+ in != RHEJ::pid::proton && in != RHEJ::pid::p_bar
+ && !RHEJ::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]);
+ }
+ result.njets += std::stoi(particles[i]);
+ }
+ else{
+ const auto pid = RHEJ::to_ParticleID(particles[i]);
+ if(!RHEJ::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;
+ }
+
+ Config to_Config(YAML::Node const & yaml){
+ try{
+ RHEJ::assert_all_options_known(yaml, get_supported_options());
+ }
+ catch(RHEJ::unknown_option const & ex){
+ throw RHEJ::unknown_option{std::string{"Unknown option '"} + ex.what() + "'"};
+ }
+
+ Config config;
+ config.process = get_process(yaml, "process");
+ set_from_yaml(config.trials, yaml, "trials");
+ 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 RHEJ::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, yaml, "unordered");
+ set_from_yaml_if_defined(config.analysis_parameters, yaml, "analysis");
+ config.scale_gen = RHEJ::ScaleGenerator{RHEJ::to_ScaleConfig(yaml)};
+ set_from_yaml_if_defined(config.output, yaml, "event output");
+ set_from_yaml_if_defined(config.RanLux_init, yaml, "RanLux init");
+ config.Higgs_coupling = RHEJ::get_Higgs_coupling(yaml, "Higgs coupling");
+ 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;
+ }
+ }
+}
diff --git a/FixedOrderGen/src/main.cc b/FixedOrderGen/src/main.cc
new file mode 100644
index 0000000..c954889
--- /dev/null
+++ b/FixedOrderGen/src/main.cc
@@ -0,0 +1,195 @@
+/**
+ * Name: main.cc
+ * Authors: Jeppe R. Andersen
+ */
+
+#include <fstream>
+#include <algorithm>
+#include <memory>
+#include <chrono>
+#include <iostream>
+#include <TROOT.h>
+#include <TFile.h>
+#include <TH1.h>
+
+#include "yaml-cpp/yaml.h"
+
+#include "config.hh"
+
+#include "LHEF/LHEF.h"
+#include "RHEJ/CombinedEventWriter.hh"
+#include "RHEJ/get_analysis.hh"
+#include "RHEJ/utility.hh"
+#include "RHEJ/PDF.hh"
+//#include "RHEJ/EventReweighter.hh"
+#include "RHEJ/stream.hh"
+#include "RHEJ/MatrixElement.hh"
+#include "RHEJ/LesHouchesWriter.hh"
+#include "PhaseSpacePoint.hh"
+
+namespace{
+ constexpr auto banner =
+ " __ ___ __ ______ __ __ \n"
+ " / / / (_)___ _/ /_ / ____/___ ___ _________ ___ __ / /__ / /______ \n"
+ " / /_/ / / __ `/ __ \\ / __/ / __ \\/ _ \\/ ___/ __ `/ / / / __ / / _ \\/ __/ ___/ \n"
+ " / __ / / /_/ / / / / / /___/ / / / __/ / / /_/ / /_/ / / /_/ / __/ /_(__ ) \n"
+ " /_/ /_/_/\\__, /_/ /_/ /_____/_/ /_/\\___/_/ \\__, /\\__, / \\____/\\___/\\__/____/ \n"
+ " ____///__/ __ ____ ///__//____/ ______ __ \n"
+ " / ____(_) _____ ____/ / / __ \\_________/ /__ _____ / ____/__ ____ ___ _________ _/ /_____ _____\n"
+ " / /_ / / |/_/ _ \\/ __ / / / / / ___/ __ / _ \\/ ___/ / / __/ _ \\/ __ \\/ _ \\/ ___/ __ `/ __/ __ \\/ ___/\n"
+ " / __/ / /> </ __/ /_/ / / /_/ / / / /_/ / __/ / / /_/ / __/ / / / __/ / / /_/ / /_/ /_/ / / \n"
+ " /_/ /_/_/|_|\\___/\\__,_/ \\____/_/ \\__,_/\\___/_/ \\____/\\___/_/ /_/\\___/_/ \\__,_/\\__/\\____/_/ \n"
+ ;
+
+ constexpr double invGeV2_to_pb = 389379292.;
+}
+
+HEJFOG::Config load_config(char const * filename){
+ try{
+ return HEJFOG::load_config(filename);
+ }
+ catch(std::exception const & exc){
+ std::cerr << "Error: " << exc.what() << '\n';
+ std::exit(EXIT_FAILURE);
+ }
+}
+
+std::unique_ptr<RHEJ::Analysis> get_analysis(
+ YAML::Node const & parameters
+){
+ try{
+ return RHEJ::get_analysis(parameters);
+ }
+ catch(std::exception const & exc){
+ std::cerr << "Failed to load analysis: " << exc.what() << '\n';
+ std::exit(EXIT_FAILURE);
+ }
+}
+
+std::string progress_bar(){
+ std::string result = "0% ";
+ for(int i = 10; i <= 100; i+= 10){
+ result += " " + std::to_string(i) + "%";
+ }
+ result += "\n|";
+ for(int i = 10; i <= 100; i+= 10){
+ result += "---------|";
+ }
+ return result + '\n';
+}
+
+int main(int argn, char** argv) {
+ if (argn < 2) {
+ std::cerr << "\n# Usage:\n.FOgen config_file\n";
+ return EXIT_FAILURE;
+ }
+ std::cout << banner;
+ fastjet::ClusterSequence::print_banner();
+ using clock = std::chrono::system_clock;
+
+ const auto start_time = clock::now();
+
+ // read configuration
+ const auto config = load_config(argv[1]);
+
+ std::unique_ptr<RHEJ::Analysis> analysis = get_analysis(
+ config.analysis_parameters
+ );
+ assert(analysis != nullptr);
+
+ RHEJ::PDF pdf{
+ config.pdf_id,
+ config.beam.particles[0], config.beam.particles[1]
+ };
+
+ // Generate a matrix element:
+ RHEJ::MatrixElement ME{
+ [&pdf](double mu){ return pdf.Halphas(mu); },
+ config.jets.def, config.jets.min_pt,
+ false,
+ config.Higgs_coupling
+ };
+
+ //TODO: fix Les Houches init block (HEPRUP)
+ LHEF::HEPRUP lhefheprup;
+ lhefheprup.IDBMUP=std::pair<long,long>(config.beam.particles[0], config.beam.particles[1]);
+ lhefheprup.EBMUP=std::make_pair(config.beam.energy, config.beam.energy);
+ lhefheprup.PDFGUP=std::make_pair(0,0);
+ lhefheprup.PDFSUP=std::make_pair(config.pdf_id,config.pdf_id);
+ lhefheprup.NPRUP=1;
+ lhefheprup.XSECUP=std::vector<double>(1.);
+ lhefheprup.XERRUP=std::vector<double>(1.);
+ lhefheprup.LPRUP=std::vector<int>{1};
+
+ RHEJ::CombinedEventWriter writer{config.output, lhefheprup};
+
+ if(config.RanLux_init){
+ HEJFOG::PhaseSpacePoint::reset_ranlux(*config.RanLux_init);
+ }
+
+ // Book root histogram for status
+ TFile hfile("GenStatus.root","RECREATE","Generator status");
+ TROOT simple("FOG","Output from HEJFOG");
+ TH1D *hmstatus;
+ hmstatus = new TH1D("mstatus","mstatus",25,-.5,24.5);
+
+
+ // Perform N trial generations
+ int iprint = 0;
+ std::cout << '\n' << progress_bar();
+ for (int trials = 0; trials< config.trials;trials++){
+ if (trials==iprint) {
+ std::cout << ".";
+ std::cout.flush();
+ iprint+=(int) config.trials/100;
+ }
+
+ // Generate phase space point
+ HEJFOG::PhaseSpacePoint psp{
+ config.process,
+ config.jets.def,config.jets.min_pt, config.jets.max_y,
+ pdf, config.beam.energy
+ };
+ hmstatus->Fill(psp.status());
+ if (psp.status()!=0) continue;
+
+ RHEJ::Event ev = config.scale_gen(
+ RHEJ::Event{
+ to_UnclusteredEvent(std::move(psp)),
+ config.jets.def, config.jets.min_pt
+ }
+ );
+ const double shat = RHEJ::shat(ev);
+ const double xa = (ev.incoming()[0].E()-ev.incoming()[0].pz())/(2.*config.beam.energy);
+ const double xb = (ev.incoming()[1].E()+ev.incoming()[1].pz())/(2.*config.beam.energy);
+
+ // evaluate matrix element on this point
+ ev.central().weight *= ME.tree(
+ ev.central().mur, ev.incoming(), ev.outgoing(), false
+ )/(shat*shat);
+ ev.central().weight *= pdf.pdfpt(0,xa,ev.central().muf, ev.incoming()[0].type);
+ ev.central().weight *= pdf.pdfpt(0,xb,ev.central().muf, ev.incoming()[1].type);
+ ev.central().weight *= invGeV2_to_pb;
+ ev.central().weight /= config.trials;
+ for(auto & var: ev.variations()){
+ var.weight *= ME.tree(
+ var.mur, ev.incoming(), ev.outgoing(), false
+ )/(shat*shat);
+ var.weight *= pdf.pdfpt(0,xa,var.muf, ev.incoming()[0].type);
+ var.weight *= pdf.pdfpt(0,xb,var.muf, ev.incoming()[1].type);
+ var.weight *= invGeV2_to_pb;
+ var.weight /= config.trials;
+ }
+ if(analysis->pass_cuts(ev)){
+ analysis->fill(ev);
+ writer.write(ev);
+ }
+
+ } // main event loop
+ std::cout << std::endl;
+ std::chrono::duration<double> run_time = (clock::now() - start_time);
+ std::cout << "\nTask Runtime: " << run_time.count() << " seconds.\n";
+
+ hfile.Write();
+ hfile.Close();
+}
diff --git a/FixedOrderGen/t/config_h_2j.yml b/FixedOrderGen/t/config_h_2j.yml
new file mode 100644
index 0000000..06c5a4c
--- /dev/null
+++ b/FixedOrderGen/t/config_h_2j.yml
@@ -0,0 +1,19 @@
+trials : 100000
+
+jets:
+ min pt: 30
+ R: 0.4
+ algorithm: antikt
+ max rapidity: 5
+
+beam:
+ energy: 6500
+ particles: [p, p]
+
+pdf: 11000
+
+process: p p => h 2j
+
+unordered: false
+
+scales: 125
diff --git a/FixedOrderGen/t/h_2j.cc b/FixedOrderGen/t/h_2j.cc
new file mode 100644
index 0000000..c2993fb
--- /dev/null
+++ b/FixedOrderGen/t/h_2j.cc
@@ -0,0 +1,79 @@
+#ifdef NDEBUG
+#undef NDEBUG
+#endif
+
+#include <algorithm>
+#include <cmath>
+#include <cassert>
+#include <iostream>
+
+#include "config.hh"
+#include "PhaseSpacePoint.hh"
+
+#include "RHEJ/PDF.hh"
+#include "RHEJ/MatrixElement.hh"
+
+using namespace HEJFOG;
+
+int main(){
+ constexpr double invGeV2_to_pb = 389379292.;
+ constexpr double xs_ref = 2.0304; //calculated with "old" HEJ svn r3364
+
+ const auto config = load_config("config_h_2j.yml");
+
+ RHEJ::PDF pdf{
+ config.pdf_id,
+ config.beam.particles[0], config.beam.particles[1]
+ };
+
+ RHEJ::MatrixElement ME{
+ [&pdf](double mu){ return pdf.Halphas(mu); },
+ config.jets.def, config.jets.min_pt,
+ false,
+ config.Higgs_coupling
+ };
+
+ double xs = 0., xs_err = 0.;
+ for (int trials = 0; trials < config.trials; ++trials){
+ HEJFOG::PhaseSpacePoint psp{
+ config.process,
+ config.jets.def,config.jets.min_pt, config.jets.max_y,
+ pdf, config.beam.energy
+ };
+ if (psp.status()!=0) continue;
+
+ RHEJ::Event ev = config.scale_gen(
+ RHEJ::Event{
+ to_UnclusteredEvent(std::move(psp)),
+ config.jets.def, config.jets.min_pt
+ }
+ );
+ const double shat = RHEJ::shat(ev);
+ const double xa = (ev.incoming()[0].E()-ev.incoming()[0].pz())/(2.*config.beam.energy);
+ const double xb = (ev.incoming()[1].E()+ev.incoming()[1].pz())/(2.*config.beam.energy);
+
+ // evaluate matrix element on this point
+ ev.central().weight *= ME.tree(
+ ev.central().mur, ev.incoming(), ev.outgoing(), false
+ )/(shat*shat);
+ ev.central().weight *= pdf.pdfpt(0,xa,ev.central().muf, ev.incoming()[0].type);
+ ev.central().weight *= pdf.pdfpt(0,xb,ev.central().muf, ev.incoming()[1].type);
+ ev.central().weight *= invGeV2_to_pb;
+ ev.central().weight /= config.trials;
+
+ const auto the_Higgs = std::find_if(
+ begin(ev.outgoing()), end(ev.outgoing()),
+ [](RHEJ::Sparticle const & p){ return p.type == RHEJ::ParticleID::h; }
+ );
+ assert(the_Higgs != end(ev.outgoing()));
+ if(std::abs(the_Higgs->rapidity()) > 5.) continue;
+
+ xs += ev.central().weight;
+ xs_err += ev.central().weight*ev.central().weight;
+ }
+ xs_err = std::sqrt(xs_err);
+ std::cout << xs_ref << " ~ " << xs << " +- " << xs_err << '\n';
+
+ assert(std::abs(xs - xs_ref) < 3*xs_err);
+ assert(xs_err < 0.01*xs);
+}
diff --git a/FixedOrderGen/t/h_3j.cc b/FixedOrderGen/t/h_3j.cc
new file mode 100644
index 0000000..e169fbf
--- /dev/null
+++ b/FixedOrderGen/t/h_3j.cc
@@ -0,0 +1,80 @@
+#ifdef NDEBUG
+#undef NDEBUG
+#endif
+
+#include <algorithm>
+#include <cmath>
+#include <cassert>
+#include <iostream>
+
+#include "config.hh"
+#include "PhaseSpacePoint.hh"
+
+#include "RHEJ/PDF.hh"
+#include "RHEJ/MatrixElement.hh"
+
+using namespace HEJFOG;
+
+int main(){
+ constexpr double invGeV2_to_pb = 389379292.;
+ constexpr double xs_ref = 1.0888; //calculated with "old" HEJ svn r3364
+
+ auto config = load_config("config_h_2j.yml");
+ config.process.njets = 3;
+
+ RHEJ::PDF pdf{
+ config.pdf_id,
+ config.beam.particles[0], config.beam.particles[1]
+ };
+
+ RHEJ::MatrixElement ME{
+ [&pdf](double mu){ return pdf.Halphas(mu); },
+ config.jets.def, config.jets.min_pt,
+ false,
+ config.Higgs_coupling
+ };
+
+ double xs = 0., xs_err = 0.;
+ for (int trials = 0; trials < config.trials; ++trials){
+ HEJFOG::PhaseSpacePoint psp{
+ config.process,
+ config.jets.def,config.jets.min_pt, config.jets.max_y,
+ pdf, config.beam.energy
+ };
+ if (psp.status()!=0) continue;
+
+ RHEJ::Event ev = config.scale_gen(
+ RHEJ::Event{
+ to_UnclusteredEvent(std::move(psp)),
+ config.jets.def, config.jets.min_pt
+ }
+ );
+ const double shat = RHEJ::shat(ev);
+ const double xa = (ev.incoming()[0].E()-ev.incoming()[0].pz())/(2.*config.beam.energy);
+ const double xb = (ev.incoming()[1].E()+ev.incoming()[1].pz())/(2.*config.beam.energy);
+
+ // evaluate matrix element on this point
+ ev.central().weight *= ME.tree(
+ ev.central().mur, ev.incoming(), ev.outgoing(), false
+ )/(shat*shat);
+ ev.central().weight *= pdf.pdfpt(0,xa,ev.central().muf, ev.incoming()[0].type);
+ ev.central().weight *= pdf.pdfpt(0,xb,ev.central().muf, ev.incoming()[1].type);
+ ev.central().weight *= invGeV2_to_pb;
+ ev.central().weight /= config.trials;
+
+ const auto the_Higgs = std::find_if(
+ begin(ev.outgoing()), end(ev.outgoing()),
+ [](RHEJ::Sparticle const & p){ return p.type == RHEJ::ParticleID::h; }
+ );
+ assert(the_Higgs != end(ev.outgoing()));
+ if(std::abs(the_Higgs->rapidity()) > 5.) continue;
+
+ xs += ev.central().weight;
+ xs_err += ev.central().weight*ev.central().weight;
+ }
+ xs_err = std::sqrt(xs_err);
+ std::cout << xs_ref << " ~ " << xs << " +- " << xs_err << '\n';
+
+ assert(std::abs(xs - xs_ref) < 3*xs_err);
+ assert(xs_err < 0.02*xs);
+}
diff --git a/analysis-plugins/CMakeLists.txt b/analysis-plugins/CMakeLists.txt
index 4af41df..3c0db33 100644
--- a/analysis-plugins/CMakeLists.txt
+++ b/analysis-plugins/CMakeLists.txt
@@ -1,14 +1,16 @@
cmake_minimum_required(VERSION 2.8 FATAL_ERROR)
set(CMAKE_EXPORT_COMPILE_COMMANDS ON)
# by default, don't expose symbols defined in plugins
set(CMAKE_CXX_VISIBILITY_PRESET hidden)
+set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} "${CMAKE_SOURCE_DIR}/../cmake/Modules/")
+
find_package(ROOT MODULE COMPONENTS Hist Tree MathCore)
include_directories(SYSTEM ${ROOT_INCLUDE_DIR})
link_directories(${ROOT_LIBRARY_DIR})
set(libraries ${ROOT_LIBRARIES} ${FASTJET_LIBRARIES} yaml-cpp rhej)
include_directories(${CMAKE_CURRENT_SOURCE_DIR}/include)
add_subdirectory(src)
diff --git a/analysis-plugins/include/StandardAnalysis_decl.hh b/analysis-plugins/include/StandardAnalysis_decl.hh
index 9b4bf5f..9a4b652 100644
--- a/analysis-plugins/include/StandardAnalysis_decl.hh
+++ b/analysis-plugins/include/StandardAnalysis_decl.hh
@@ -1,143 +1,146 @@
#pragma once
#include <array>
#include <map>
#include "fastjet/JetDefinition.hh"
#include "fastjet/PseudoJet.hh"
#include "fastjet/ClusterSequence.hh"
#include "RHEJ/Analysis.hh"
#include "RHEJ/Event.hh"
#include "Histogram.hh"
class TFile;
namespace YAML{
class Node;
}
namespace RHEJ_analyses{
using RHEJ::Event;
namespace detail{
struct AbsWtLess{
bool operator()(Event const & a, Event const & b){
return std::abs(a.central().weight) < std::abs(b.central().weight);
}
};
}
class StandardAnalysis: public RHEJ::DynamicAnalysis{
public:
using histogram_type = Histogram;
// macro definition to avoid repetition of histogram names
#define RHEJ_HISTOGRAMS(F) F(total), F(np), F(njets),F(pt1), F(pt2), F(pth1), F(pth2), F(qt), F(ya), F(yb), F(ydif), F(yEW), F(ptEW), F(wt), F(wtwt)
#define RHEJ_AS_ENUM(VAR) VAR
#define RHEJ_AS_STRING(VAR) #VAR
enum Id{
RHEJ_HISTOGRAMS(RHEJ_AS_ENUM),
first_histogram = np,
last_histogram = wtwt,
first_histpack = np,
last_histpack = ptEW
};
static constexpr auto hist_names = RHEJ::make_array(
RHEJ_HISTOGRAMS(RHEJ_AS_STRING)
);
#undef RHEJ_HISTOGRAMS
#undef RHEJ_AS_ENUM
#undef RHEJ_AS_STRING
StandardAnalysis(YAML::Node const &);
explicit StandardAnalysis(std::string const & filename);
StandardAnalysis(StandardAnalysis const & other);
StandardAnalysis(StandardAnalysis && other);
StandardAnalysis & operator=(StandardAnalysis const & other);
StandardAnalysis & operator=(StandardAnalysis && other);
void fill(Event const & e) override;
+ bool pass_cuts(Event const &) override{
+ return true;
+ };
std::vector<std::string> histogram_names() const override;
std::vector<Parameter> parameters() const override;
void draw_histogram(
std::string const & name, std::string const & filename
) const override;
bool update_parameter(std::string const & name, double value) override;
~StandardAnalysis() override;
static std::unique_ptr<Analysis> create(YAML::Node const &);
static std::unique_ptr<DynamicAnalysis> create(
std::string const & filename
);
private:
void store();
void restore(std::string const & filename);
void write_settings(TFile & out);
void write_histograms(TFile & out);
void write_events(TFile & out);
void read_settings(TFile & in);
void read_histograms(TFile & in);
void read_events(TFile & in);
struct hist_pack{
histogram_type h, hwt, hN;
hist_pack() = default;
hist_pack(
const char * title, const std::string & id,
int numbins, double min, double max
);
void fill(double val, double wt);
void write() const;
int nbins() const;
};
struct KinInfo{
std::vector<RHEJ::Sparticle> out;
std::vector<fastjet::PseudoJet> jets_y;
std::vector<fastjet::PseudoJet> jets_pt;
fastjet::PseudoJet q;
};
struct EventHists{
EventHists();
explicit EventHists(std::string const & id);
void fill(KinInfo const & i, double weight);
void write() const;
std::array<hist_pack, last_histpack + 1> hist_;
histogram_type hist_wt_, hist_wtwt_;
};
void restore_histograms(
std::map<std::string, histogram_type> && hists
);
void update_hist_by_name();
void add_to_hist_by_name(EventHists const & hist);
void add_to_hist_by_name(hist_pack const & pack);
void reset();
EventHists central_;
std::vector<EventHists> hist_;
double log_wt_cut_;
std::string output_;
std::set<Event, detail::AbsWtLess> bad_events_;
using Hist_ref = std::reference_wrapper<const histogram_type>;
std::map<std::string, Hist_ref> hist_by_name_;
EventHists orig_central_;
std::vector<EventHists> orig_hist_;
double orig_log_wt_cut_;
};
}
diff --git a/analysis-plugins/src/VBF.cc b/analysis-plugins/src/VBF.cc
index 95e869f..3d1d045 100644
--- a/analysis-plugins/src/VBF.cc
+++ b/analysis-plugins/src/VBF.cc
@@ -1,93 +1,97 @@
#include "StandardAnalysis.hh"
#include <iostream>
namespace{
using namespace RHEJ_analyses;
class VBF: RHEJ::DynamicAnalysis{
public:
VBF(YAML::Node const & parameters):
analysis_{parameters},
min_dy12_{parameters["min dy12"].as<double>()},
min_m12_{parameters["min m12"].as<double>()}
{}
VBF(std::string const & filename):
analysis_{filename},
/**
* The following lines are only safe if "filename" was actually
* generated using this analysis.
*
* In this case all stored event already fulfil the VBF cuts
* and we can ignore them
*/
min_dy12_{0.},
min_m12_{0.}
{}
void fill(Event const & ev) override{
+ analysis_.fill(ev);
+ }
+
+ bool pass_cuts(Event const & ev) override{
const auto jets_pt = sorted_by_pt(ev.jets());
// tight cuts a la arXiv:1311.6738
- if(jets_pt.size() < 2) return;
+ if(jets_pt.size() < 2) return false;
const double y0 = jets_pt[0].rapidity();
const double y1 = jets_pt[1].rapidity();
- if(y0*y1 >= 0) return;
- if(std::abs(y0 - y1) <= min_dy12_) return;
- if((jets_pt[0] + jets_pt[1]).m2() <= min_m12_*min_m12_) return;
- analysis_.fill(ev);
+ if(y0*y1 >= 0) return false;
+ if(std::abs(y0 - y1) <= min_dy12_) return false;
+ if((jets_pt[0] + jets_pt[1]).m2() <= min_m12_*min_m12_) return false;
+ return true;
}
std::vector<std::string> histogram_names() const override{
return analysis_.histogram_names();
}
std::vector<Parameter> parameters() const override{
return analysis_.parameters();
}
void draw_histogram(
std::string const & name, std::string const & filename
) const override{
analysis_.draw_histogram(name, filename);
}
bool update_parameter(std::string const & name, double value) override{
return analysis_.update_parameter(name, value);
}
static std::unique_ptr<RHEJ::Analysis> create(
YAML::Node const & parameters
){
return std::unique_ptr<RHEJ::Analysis>{new VBF{parameters}};
}
static std::unique_ptr<DynamicAnalysis> create(
std::string const & filename
){
return std::unique_ptr<RHEJ::DynamicAnalysis>{new VBF{filename}};
}
private:
StandardAnalysis analysis_;
double min_dy12_;
double min_m12_;
};
}
extern "C"
__attribute__((visibility("default")))
std::unique_ptr<RHEJ::Analysis> make_analysis(
YAML::Node const & parameters
){
return VBF::create(parameters);
}
extern "C"
__attribute__((visibility("default")))
std::unique_ptr<RHEJ::Analysis> load_analysis(
std::string const & filename
){
return VBF::create(filename);
}
diff --git a/doc/.gitignore b/doc/.gitignore
new file mode 100644
index 0000000..fe056cc
--- /dev/null
+++ b/doc/.gitignore
@@ -0,0 +1,2 @@
+html/
+latex/
\ No newline at end of file
diff --git a/doc/doxygen/Doxyfile b/doc/doxygen/Doxyfile
new file mode 100644
index 0000000..b29f0f9
--- /dev/null
+++ b/doc/doxygen/Doxyfile
@@ -0,0 +1,2494 @@
+# Doxyfile 1.8.13
+
+# This file describes the settings to be used by the documentation system
+# doxygen (www.doxygen.org) for a project.
+#
+# All text after a double hash (##) is considered a comment and is placed in
+# front of the TAG it is preceding.
+#
+# All text after a single hash (#) is considered a comment and will be ignored.
+# The format is:
+# TAG = value [value, ...]
+# For lists, items can also be appended using:
+# TAG += value [value, ...]
+# Values that contain spaces should be placed between quotes (\" \").
+
+#---------------------------------------------------------------------------
+# Project related configuration options
+#---------------------------------------------------------------------------
+
+# This tag specifies the encoding used for all characters in the config file
+# that follow. The default is UTF-8 which is also the encoding used for all text
+# before the first occurrence of this tag. Doxygen uses libiconv (or the iconv
+# built into libc) for the transcoding. See http://www.gnu.org/software/libiconv
+# for the list of possible encodings.
+# The default value is: UTF-8.
+
+DOXYFILE_ENCODING = UTF-8
+
+# The PROJECT_NAME tag is a single word (or a sequence of words surrounded by
+# double-quotes, unless you are using Doxywizard) that should identify the
+# project for which the documentation is generated. This name is used in the
+# title of most generated pages and in a few other places.
+# The default value is: My Project.
+
+PROJECT_NAME = "My Project"
+
+# The PROJECT_NUMBER tag can be used to enter a project or revision number. This
+# could be handy for archiving the generated documentation or if some version
+# control system is used.
+
+PROJECT_NUMBER =
+
+# Using the PROJECT_BRIEF tag one can provide an optional one line description
+# for a project that appears at the top of each page and should give viewer a
+# quick idea about the purpose of the project. Keep the description short.
+
+PROJECT_BRIEF =
+
+# With the PROJECT_LOGO tag one can specify a logo or an icon that is included
+# in the documentation. The maximum height of the logo should not exceed 55
+# pixels and the maximum width should not exceed 200 pixels. Doxygen will copy
+# the logo to the output directory.
+
+PROJECT_LOGO =
+
+# The OUTPUT_DIRECTORY tag is used to specify the (relative or absolute) path
+# into which the generated documentation will be written. If a relative path is
+# entered, it will be relative to the location where doxygen was started. If
+# left blank the current directory will be used.
+
+OUTPUT_DIRECTORY =
+
+# If the CREATE_SUBDIRS tag is set to YES then doxygen will create 4096 sub-
+# directories (in 2 levels) under the output directory of each output format and
+# will distribute the generated files over these directories. Enabling this
+# option can be useful when feeding doxygen a huge amount of source files, where
+# putting all generated files in the same directory would otherwise causes
+# performance problems for the file system.
+# The default value is: NO.
+
+CREATE_SUBDIRS = NO
+
+# If the ALLOW_UNICODE_NAMES tag is set to YES, doxygen will allow non-ASCII
+# characters to appear in the names of generated files. If set to NO, non-ASCII
+# characters will be escaped, for example _xE3_x81_x84 will be used for Unicode
+# U+3044.
+# The default value is: NO.
+
+ALLOW_UNICODE_NAMES = NO
+
+# The OUTPUT_LANGUAGE tag is used to specify the language in which all
+# documentation generated by doxygen is written. Doxygen will use this
+# information to generate all constant output in the proper language.
+# Possible values are: Afrikaans, Arabic, Armenian, Brazilian, Catalan, Chinese,
+# Chinese-Traditional, Croatian, Czech, Danish, Dutch, English (United States),
+# Esperanto, Farsi (Persian), Finnish, French, German, Greek, Hungarian,
+# Indonesian, Italian, Japanese, Japanese-en (Japanese with English messages),
+# Korean, Korean-en (Korean with English messages), Latvian, Lithuanian,
+# Macedonian, Norwegian, Persian (Farsi), Polish, Portuguese, Romanian, Russian,
+# Serbian, Serbian-Cyrillic, Slovak, Slovene, Spanish, Swedish, Turkish,
+# Ukrainian and Vietnamese.
+# The default value is: English.
+
+OUTPUT_LANGUAGE = English
+
+# If the BRIEF_MEMBER_DESC tag is set to YES, doxygen will include brief member
+# descriptions after the members that are listed in the file and class
+# documentation (similar to Javadoc). Set to NO to disable this.
+# The default value is: YES.
+
+BRIEF_MEMBER_DESC = YES
+
+# If the REPEAT_BRIEF tag is set to YES, doxygen will prepend the brief
+# description of a member or function before the detailed description
+#
+# Note: If both HIDE_UNDOC_MEMBERS and BRIEF_MEMBER_DESC are set to NO, the
+# brief descriptions will be completely suppressed.
+# The default value is: YES.
+
+REPEAT_BRIEF = YES
+
+# This tag implements a quasi-intelligent brief description abbreviator that is
+# used to form the text in various listings. Each string in this list, if found
+# as the leading text of the brief description, will be stripped from the text
+# and the result, after processing the whole list, is used as the annotated
+# text. Otherwise, the brief description is used as-is. If left blank, the
+# following values are used ($name is automatically replaced with the name of
+# the entity):The $name class, The $name widget, The $name file, is, provides,
+# specifies, contains, represents, a, an and the.
+
+ABBREVIATE_BRIEF = "The $name class" \
+ "The $name widget" \
+ "The $name file" \
+ is \
+ provides \
+ specifies \
+ contains \
+ represents \
+ a \
+ an \
+ the
+
+# If the ALWAYS_DETAILED_SEC and REPEAT_BRIEF tags are both set to YES then
+# doxygen will generate a detailed section even if there is only a brief
+# description.
+# The default value is: NO.
+
+ALWAYS_DETAILED_SEC = NO
+
+# If the INLINE_INHERITED_MEMB tag is set to YES, doxygen will show all
+# inherited members of a class in the documentation of that class as if those
+# members were ordinary class members. Constructors, destructors and assignment
+# operators of the base classes will not be shown.
+# The default value is: NO.
+
+INLINE_INHERITED_MEMB = NO
+
+# If the FULL_PATH_NAMES tag is set to YES, doxygen will prepend the full path
+# before files name in the file list and in the header files. If set to NO the
+# shortest path that makes the file name unique will be used
+# The default value is: YES.
+
+FULL_PATH_NAMES = YES
+
+# The STRIP_FROM_PATH tag can be used to strip a user-defined part of the path.
+# Stripping is only done if one of the specified strings matches the left-hand
+# part of the path. The tag can be used to show relative paths in the file list.
+# If left blank the directory from which doxygen is run is used as the path to
+# strip.
+#
+# Note that you can specify absolute paths here, but also relative paths, which
+# will be relative from the directory where doxygen is started.
+# This tag requires that the tag FULL_PATH_NAMES is set to YES.
+
+STRIP_FROM_PATH =
+
+# The STRIP_FROM_INC_PATH tag can be used to strip a user-defined part of the
+# path mentioned in the documentation of a class, which tells the reader which
+# header file to include in order to use a class. If left blank only the name of
+# the header file containing the class definition is used. Otherwise one should
+# specify the list of include paths that are normally passed to the compiler
+# using the -I flag.
+
+STRIP_FROM_INC_PATH =
+
+# If the SHORT_NAMES tag is set to YES, doxygen will generate much shorter (but
+# less readable) file names. This can be useful is your file systems doesn't
+# support long names like on DOS, Mac, or CD-ROM.
+# The default value is: NO.
+
+SHORT_NAMES = NO
+
+# If the JAVADOC_AUTOBRIEF tag is set to YES then doxygen will interpret the
+# first line (until the first dot) of a Javadoc-style comment as the brief
+# description. If set to NO, the Javadoc-style will behave just like regular Qt-
+# style comments (thus requiring an explicit @brief command for a brief
+# description.)
+# The default value is: NO.
+
+JAVADOC_AUTOBRIEF = NO
+
+# If the QT_AUTOBRIEF tag is set to YES then doxygen will interpret the first
+# line (until the first dot) of a Qt-style comment as the brief description. If
+# set to NO, the Qt-style will behave just like regular Qt-style comments (thus
+# requiring an explicit \brief command for a brief description.)
+# The default value is: NO.
+
+QT_AUTOBRIEF = NO
+
+# The MULTILINE_CPP_IS_BRIEF tag can be set to YES to make doxygen treat a
+# multi-line C++ special comment block (i.e. a block of //! or /// comments) as
+# a brief description. This used to be the default behavior. The new default is
+# to treat a multi-line C++ comment block as a detailed description. Set this
+# tag to YES if you prefer the old behavior instead.
+#
+# Note that setting this tag to YES also means that rational rose comments are
+# not recognized any more.
+# The default value is: NO.
+
+MULTILINE_CPP_IS_BRIEF = NO
+
+# If the INHERIT_DOCS tag is set to YES then an undocumented member inherits the
+# documentation from any documented member that it re-implements.
+# The default value is: YES.
+
+INHERIT_DOCS = YES
+
+# If the SEPARATE_MEMBER_PAGES tag is set to YES then doxygen will produce a new
+# page for each member. If set to NO, the documentation of a member will be part
+# of the file/class/namespace that contains it.
+# The default value is: NO.
+
+SEPARATE_MEMBER_PAGES = NO
+
+# The TAB_SIZE tag can be used to set the number of spaces in a tab. Doxygen
+# uses this value to replace tabs by spaces in code fragments.
+# Minimum value: 1, maximum value: 16, default value: 4.
+
+TAB_SIZE = 4
+
+# This tag can be used to specify a number of aliases that act as commands in
+# the documentation. An alias has the form:
+# name=value
+# For example adding
+# "sideeffect=@par Side Effects:\n"
+# will allow you to put the command \sideeffect (or @sideeffect) in the
+# documentation, which will result in a user-defined paragraph with heading
+# "Side Effects:". You can put \n's in the value part of an alias to insert
+# newlines.
+
+ALIASES =
+
+# This tag can be used to specify a number of word-keyword mappings (TCL only).
+# A mapping has the form "name=value". For example adding "class=itcl::class"
+# will allow you to use the command class in the itcl::class meaning.
+
+TCL_SUBST =
+
+# Set the OPTIMIZE_OUTPUT_FOR_C tag to YES if your project consists of C sources
+# only. Doxygen will then generate output that is more tailored for C. For
+# instance, some of the names that are used will be different. The list of all
+# members will be omitted, etc.
+# The default value is: NO.
+
+OPTIMIZE_OUTPUT_FOR_C = NO
+
+# Set the OPTIMIZE_OUTPUT_JAVA tag to YES if your project consists of Java or
+# Python sources only. Doxygen will then generate output that is more tailored
+# for that language. For instance, namespaces will be presented as packages,
+# qualified scopes will look different, etc.
+# The default value is: NO.
+
+OPTIMIZE_OUTPUT_JAVA = NO
+
+# Set the OPTIMIZE_FOR_FORTRAN tag to YES if your project consists of Fortran
+# sources. Doxygen will then generate output that is tailored for Fortran.
+# The default value is: NO.
+
+OPTIMIZE_FOR_FORTRAN = NO
+
+# Set the OPTIMIZE_OUTPUT_VHDL tag to YES if your project consists of VHDL
+# sources. Doxygen will then generate output that is tailored for VHDL.
+# The default value is: NO.
+
+OPTIMIZE_OUTPUT_VHDL = NO
+
+# Doxygen selects the parser to use depending on the extension of the files it
+# parses. With this tag you can assign which parser to use for a given
+# extension. Doxygen has a built-in mapping, but you can override or extend it
+# using this tag. The format is ext=language, where ext is a file extension, and
+# language is one of the parsers supported by doxygen: IDL, Java, Javascript,
+# C#, C, C++, D, PHP, Objective-C, Python, Fortran (fixed format Fortran:
+# FortranFixed, free formatted Fortran: FortranFree, unknown formatted Fortran:
+# Fortran. In the later case the parser tries to guess whether the code is fixed
+# or free formatted code, this is the default for Fortran type files), VHDL. For
+# instance to make doxygen treat .inc files as Fortran files (default is PHP),
+# and .f files as C (default is Fortran), use: inc=Fortran f=C.
+#
+# Note: For files without extension you can use no_extension as a placeholder.
+#
+# Note that for custom extensions you also need to set FILE_PATTERNS otherwise
+# the files are not read by doxygen.
+
+EXTENSION_MAPPING =
+
+# If the MARKDOWN_SUPPORT tag is enabled then doxygen pre-processes all comments
+# according to the Markdown format, which allows for more readable
+# documentation. See http://daringfireball.net/projects/markdown/ for details.
+# The output of markdown processing is further processed by doxygen, so you can
+# mix doxygen, HTML, and XML commands with Markdown formatting. Disable only in
+# case of backward compatibilities issues.
+# The default value is: YES.
+
+MARKDOWN_SUPPORT = YES
+
+# When the TOC_INCLUDE_HEADINGS tag is set to a non-zero value, all headings up
+# to that level are automatically included in the table of contents, even if
+# they do not have an id attribute.
+# Note: This feature currently applies only to Markdown headings.
+# Minimum value: 0, maximum value: 99, default value: 0.
+# This tag requires that the tag MARKDOWN_SUPPORT is set to YES.
+
+TOC_INCLUDE_HEADINGS = 0
+
+# When enabled doxygen tries to link words that correspond to documented
+# classes, or namespaces to their corresponding documentation. Such a link can
+# be prevented in individual cases by putting a % sign in front of the word or
+# globally by setting AUTOLINK_SUPPORT to NO.
+# The default value is: YES.
+
+AUTOLINK_SUPPORT = YES
+
+# If you use STL classes (i.e. std::string, std::vector, etc.) but do not want
+# to include (a tag file for) the STL sources as input, then you should set this
+# tag to YES in order to let doxygen match functions declarations and
+# definitions whose arguments contain STL classes (e.g. func(std::string);
+# versus func(std::string) {}). This also make the inheritance and collaboration
+# diagrams that involve STL classes more complete and accurate.
+# The default value is: NO.
+
+BUILTIN_STL_SUPPORT = NO
+
+# If you use Microsoft's C++/CLI language, you should set this option to YES to
+# enable parsing support.
+# The default value is: NO.
+
+CPP_CLI_SUPPORT = NO
+
+# Set the SIP_SUPPORT tag to YES if your project consists of sip (see:
+# http://www.riverbankcomputing.co.uk/software/sip/intro) sources only. Doxygen
+# will parse them like normal C++ but will assume all classes use public instead
+# of private inheritance when no explicit protection keyword is present.
+# The default value is: NO.
+
+SIP_SUPPORT = NO
+
+# For Microsoft's IDL there are propget and propput attributes to indicate
+# getter and setter methods for a property. Setting this option to YES will make
+# doxygen to replace the get and set methods by a property in the documentation.
+# This will only work if the methods are indeed getting or setting a simple
+# type. If this is not the case, or you want to show the methods anyway, you
+# should set this option to NO.
+# The default value is: YES.
+
+IDL_PROPERTY_SUPPORT = YES
+
+# If member grouping is used in the documentation and the DISTRIBUTE_GROUP_DOC
+# tag is set to YES then doxygen will reuse the documentation of the first
+# member in the group (if any) for the other members of the group. By default
+# all members of a group must be documented explicitly.
+# The default value is: NO.
+
+DISTRIBUTE_GROUP_DOC = NO
+
+# If one adds a struct or class to a group and this option is enabled, then also
+# any nested class or struct is added to the same group. By default this option
+# is disabled and one has to add nested compounds explicitly via \ingroup.
+# The default value is: NO.
+
+GROUP_NESTED_COMPOUNDS = NO
+
+# Set the SUBGROUPING tag to YES to allow class member groups of the same type
+# (for instance a group of public functions) to be put as a subgroup of that
+# type (e.g. under the Public Functions section). Set it to NO to prevent
+# subgrouping. Alternatively, this can be done per class using the
+# \nosubgrouping command.
+# The default value is: YES.
+
+SUBGROUPING = YES
+
+# When the INLINE_GROUPED_CLASSES tag is set to YES, classes, structs and unions
+# are shown inside the group in which they are included (e.g. using \ingroup)
+# instead of on a separate page (for HTML and Man pages) or section (for LaTeX
+# and RTF).
+#
+# Note that this feature does not work in combination with
+# SEPARATE_MEMBER_PAGES.
+# The default value is: NO.
+
+INLINE_GROUPED_CLASSES = NO
+
+# When the INLINE_SIMPLE_STRUCTS tag is set to YES, structs, classes, and unions
+# with only public data fields or simple typedef fields will be shown inline in
+# the documentation of the scope in which they are defined (i.e. file,
+# namespace, or group documentation), provided this scope is documented. If set
+# to NO, structs, classes, and unions are shown on a separate page (for HTML and
+# Man pages) or section (for LaTeX and RTF).
+# The default value is: NO.
+
+INLINE_SIMPLE_STRUCTS = NO
+
+# When TYPEDEF_HIDES_STRUCT tag is enabled, a typedef of a struct, union, or
+# enum is documented as struct, union, or enum with the name of the typedef. So
+# typedef struct TypeS {} TypeT, will appear in the documentation as a struct
+# with name TypeT. When disabled the typedef will appear as a member of a file,
+# namespace, or class. And the struct will be named TypeS. This can typically be
+# useful for C code in case the coding convention dictates that all compound
+# types are typedef'ed and only the typedef is referenced, never the tag name.
+# The default value is: NO.
+
+TYPEDEF_HIDES_STRUCT = NO
+
+# The size of the symbol lookup cache can be set using LOOKUP_CACHE_SIZE. This
+# cache is used to resolve symbols given their name and scope. Since this can be
+# an expensive process and often the same symbol appears multiple times in the
+# code, doxygen keeps a cache of pre-resolved symbols. If the cache is too small
+# doxygen will become slower. If the cache is too large, memory is wasted. The
+# cache size is given by this formula: 2^(16+LOOKUP_CACHE_SIZE). The valid range
+# is 0..9, the default is 0, corresponding to a cache size of 2^16=65536
+# symbols. At the end of a run doxygen will report the cache usage and suggest
+# the optimal cache size from a speed point of view.
+# Minimum value: 0, maximum value: 9, default value: 0.
+
+LOOKUP_CACHE_SIZE = 0
+
+#---------------------------------------------------------------------------
+# Build related configuration options
+#---------------------------------------------------------------------------
+
+# If the EXTRACT_ALL tag is set to YES, doxygen will assume all entities in
+# documentation are documented, even if no documentation was available. Private
+# class members and static file members will be hidden unless the
+# EXTRACT_PRIVATE respectively EXTRACT_STATIC tags are set to YES.
+# Note: This will also disable the warnings about undocumented members that are
+# normally produced when WARNINGS is set to YES.
+# The default value is: NO.
+
+EXTRACT_ALL = NO
+
+# If the EXTRACT_PRIVATE tag is set to YES, all private members of a class will
+# be included in the documentation.
+# The default value is: NO.
+
+EXTRACT_PRIVATE = NO
+
+# If the EXTRACT_PACKAGE tag is set to YES, all members with package or internal
+# scope will be included in the documentation.
+# The default value is: NO.
+
+EXTRACT_PACKAGE = NO
+
+# If the EXTRACT_STATIC tag is set to YES, all static members of a file will be
+# included in the documentation.
+# The default value is: NO.
+
+EXTRACT_STATIC = NO
+
+# If the EXTRACT_LOCAL_CLASSES tag is set to YES, classes (and structs) defined
+# locally in source files will be included in the documentation. If set to NO,
+# only classes defined in header files are included. Does not have any effect
+# for Java sources.
+# The default value is: YES.
+
+EXTRACT_LOCAL_CLASSES = YES
+
+# This flag is only useful for Objective-C code. If set to YES, local methods,
+# which are defined in the implementation section but not in the interface are
+# included in the documentation. If set to NO, only methods in the interface are
+# included.
+# The default value is: NO.
+
+EXTRACT_LOCAL_METHODS = NO
+
+# If this flag is set to YES, the members of anonymous namespaces will be
+# extracted and appear in the documentation as a namespace called
+# 'anonymous_namespace{file}', where file will be replaced with the base name of
+# the file that contains the anonymous namespace. By default anonymous namespace
+# are hidden.
+# The default value is: NO.
+
+EXTRACT_ANON_NSPACES = NO
+
+# If the HIDE_UNDOC_MEMBERS tag is set to YES, doxygen will hide all
+# undocumented members inside documented classes or files. If set to NO these
+# members will be included in the various overviews, but no documentation
+# section is generated. This option has no effect if EXTRACT_ALL is enabled.
+# The default value is: NO.
+
+HIDE_UNDOC_MEMBERS = NO
+
+# If the HIDE_UNDOC_CLASSES tag is set to YES, doxygen will hide all
+# undocumented classes that are normally visible in the class hierarchy. If set
+# to NO, these classes will be included in the various overviews. This option
+# has no effect if EXTRACT_ALL is enabled.
+# The default value is: NO.
+
+HIDE_UNDOC_CLASSES = NO
+
+# If the HIDE_FRIEND_COMPOUNDS tag is set to YES, doxygen will hide all friend
+# (class|struct|union) declarations. If set to NO, these declarations will be
+# included in the documentation.
+# The default value is: NO.
+
+HIDE_FRIEND_COMPOUNDS = NO
+
+# If the HIDE_IN_BODY_DOCS tag is set to YES, doxygen will hide any
+# documentation blocks found inside the body of a function. If set to NO, these
+# blocks will be appended to the function's detailed documentation block.
+# The default value is: NO.
+
+HIDE_IN_BODY_DOCS = NO
+
+# The INTERNAL_DOCS tag determines if documentation that is typed after a
+# \internal command is included. If the tag is set to NO then the documentation
+# will be excluded. Set it to YES to include the internal documentation.
+# The default value is: NO.
+
+INTERNAL_DOCS = NO
+
+# If the CASE_SENSE_NAMES tag is set to NO then doxygen will only generate file
+# names in lower-case letters. If set to YES, upper-case letters are also
+# allowed. This is useful if you have classes or files whose names only differ
+# in case and if your file system supports case sensitive file names. Windows
+# and Mac users are advised to set this option to NO.
+# The default value is: system dependent.
+
+CASE_SENSE_NAMES = YES
+
+# If the HIDE_SCOPE_NAMES tag is set to NO then doxygen will show members with
+# their full class and namespace scopes in the documentation. If set to YES, the
+# scope will be hidden.
+# The default value is: NO.
+
+HIDE_SCOPE_NAMES = NO
+
+# If the HIDE_COMPOUND_REFERENCE tag is set to NO (default) then doxygen will
+# append additional text to a page's title, such as Class Reference. If set to
+# YES the compound reference will be hidden.
+# The default value is: NO.
+
+HIDE_COMPOUND_REFERENCE= NO
+
+# If the SHOW_INCLUDE_FILES tag is set to YES then doxygen will put a list of
+# the files that are included by a file in the documentation of that file.
+# The default value is: YES.
+
+SHOW_INCLUDE_FILES = YES
+
+# If the SHOW_GROUPED_MEMB_INC tag is set to YES then Doxygen will add for each
+# grouped member an include statement to the documentation, telling the reader
+# which file to include in order to use the member.
+# The default value is: NO.
+
+SHOW_GROUPED_MEMB_INC = NO
+
+# If the FORCE_LOCAL_INCLUDES tag is set to YES then doxygen will list include
+# files with double quotes in the documentation rather than with sharp brackets.
+# The default value is: NO.
+
+FORCE_LOCAL_INCLUDES = NO
+
+# If the INLINE_INFO tag is set to YES then a tag [inline] is inserted in the
+# documentation for inline members.
+# The default value is: YES.
+
+INLINE_INFO = YES
+
+# If the SORT_MEMBER_DOCS tag is set to YES then doxygen will sort the
+# (detailed) documentation of file and class members alphabetically by member
+# name. If set to NO, the members will appear in declaration order.
+# The default value is: YES.
+
+SORT_MEMBER_DOCS = YES
+
+# If the SORT_BRIEF_DOCS tag is set to YES then doxygen will sort the brief
+# descriptions of file, namespace and class members alphabetically by member
+# name. If set to NO, the members will appear in declaration order. Note that
+# this will also influence the order of the classes in the class list.
+# The default value is: NO.
+
+SORT_BRIEF_DOCS = NO
+
+# If the SORT_MEMBERS_CTORS_1ST tag is set to YES then doxygen will sort the
+# (brief and detailed) documentation of class members so that constructors and
+# destructors are listed first. If set to NO the constructors will appear in the
+# respective orders defined by SORT_BRIEF_DOCS and SORT_MEMBER_DOCS.
+# Note: If SORT_BRIEF_DOCS is set to NO this option is ignored for sorting brief
+# member documentation.
+# Note: If SORT_MEMBER_DOCS is set to NO this option is ignored for sorting
+# detailed member documentation.
+# The default value is: NO.
+
+SORT_MEMBERS_CTORS_1ST = NO
+
+# If the SORT_GROUP_NAMES tag is set to YES then doxygen will sort the hierarchy
+# of group names into alphabetical order. If set to NO the group names will
+# appear in their defined order.
+# The default value is: NO.
+
+SORT_GROUP_NAMES = NO
+
+# If the SORT_BY_SCOPE_NAME tag is set to YES, the class list will be sorted by
+# fully-qualified names, including namespaces. If set to NO, the class list will
+# be sorted only by class name, not including the namespace part.
+# Note: This option is not very useful if HIDE_SCOPE_NAMES is set to YES.
+# Note: This option applies only to the class list, not to the alphabetical
+# list.
+# The default value is: NO.
+
+SORT_BY_SCOPE_NAME = NO
+
+# If the STRICT_PROTO_MATCHING option is enabled and doxygen fails to do proper
+# type resolution of all parameters of a function it will reject a match between
+# the prototype and the implementation of a member function even if there is
+# only one candidate or it is obvious which candidate to choose by doing a
+# simple string match. By disabling STRICT_PROTO_MATCHING doxygen will still
+# accept a match between prototype and implementation in such cases.
+# The default value is: NO.
+
+STRICT_PROTO_MATCHING = NO
+
+# The GENERATE_TODOLIST tag can be used to enable (YES) or disable (NO) the todo
+# list. This list is created by putting \todo commands in the documentation.
+# The default value is: YES.
+
+GENERATE_TODOLIST = YES
+
+# The GENERATE_TESTLIST tag can be used to enable (YES) or disable (NO) the test
+# list. This list is created by putting \test commands in the documentation.
+# The default value is: YES.
+
+GENERATE_TESTLIST = YES
+
+# The GENERATE_BUGLIST tag can be used to enable (YES) or disable (NO) the bug
+# list. This list is created by putting \bug commands in the documentation.
+# The default value is: YES.
+
+GENERATE_BUGLIST = YES
+
+# The GENERATE_DEPRECATEDLIST tag can be used to enable (YES) or disable (NO)
+# the deprecated list. This list is created by putting \deprecated commands in
+# the documentation.
+# The default value is: YES.
+
+GENERATE_DEPRECATEDLIST= YES
+
+# The ENABLED_SECTIONS tag can be used to enable conditional documentation
+# sections, marked by \if <section_label> ... \endif and \cond <section_label>
+# ... \endcond blocks.
+
+ENABLED_SECTIONS =
+
+# The MAX_INITIALIZER_LINES tag determines the maximum number of lines that the
+# initial value of a variable or macro / define can have for it to appear in the
+# documentation. If the initializer consists of more lines than specified here
+# it will be hidden. Use a value of 0 to hide initializers completely. The
+# appearance of the value of individual variables and macros / defines can be
+# controlled using \showinitializer or \hideinitializer command in the
+# documentation regardless of this setting.
+# Minimum value: 0, maximum value: 10000, default value: 30.
+
+MAX_INITIALIZER_LINES = 30
+
+# Set the SHOW_USED_FILES tag to NO to disable the list of files generated at
+# the bottom of the documentation of classes and structs. If set to YES, the
+# list will mention the files that were used to generate the documentation.
+# The default value is: YES.
+
+SHOW_USED_FILES = YES
+
+# Set the SHOW_FILES tag to NO to disable the generation of the Files page. This
+# will remove the Files entry from the Quick Index and from the Folder Tree View
+# (if specified).
+# The default value is: YES.
+
+SHOW_FILES = YES
+
+# Set the SHOW_NAMESPACES tag to NO to disable the generation of the Namespaces
+# page. This will remove the Namespaces entry from the Quick Index and from the
+# Folder Tree View (if specified).
+# The default value is: YES.
+
+SHOW_NAMESPACES = YES
+
+# The FILE_VERSION_FILTER tag can be used to specify a program or script that
+# doxygen should invoke to get the current version for each file (typically from
+# the version control system). Doxygen will invoke the program by executing (via
+# popen()) the command command input-file, where command is the value of the
+# FILE_VERSION_FILTER tag, and input-file is the name of an input file provided
+# by doxygen. Whatever the program writes to standard output is used as the file
+# version. For an example see the documentation.
+
+FILE_VERSION_FILTER =
+
+# The LAYOUT_FILE tag can be used to specify a layout file which will be parsed
+# by doxygen. The layout file controls the global structure of the generated
+# output files in an output format independent way. To create the layout file
+# that represents doxygen's defaults, run doxygen with the -l option. You can
+# optionally specify a file name after the option, if omitted DoxygenLayout.xml
+# will be used as the name of the layout file.
+#
+# Note that if you run doxygen from a directory containing a file called
+# DoxygenLayout.xml, doxygen will parse it automatically even if the LAYOUT_FILE
+# tag is left empty.
+
+LAYOUT_FILE =
+
+# The CITE_BIB_FILES tag can be used to specify one or more bib files containing
+# the reference definitions. This must be a list of .bib files. The .bib
+# extension is automatically appended if omitted. This requires the bibtex tool
+# to be installed. See also http://en.wikipedia.org/wiki/BibTeX for more info.
+# For LaTeX the style of the bibliography can be controlled using
+# LATEX_BIB_STYLE. To use this feature you need bibtex and perl available in the
+# search path. See also \cite for info how to create references.
+
+CITE_BIB_FILES =
+
+#---------------------------------------------------------------------------
+# Configuration options related to warning and progress messages
+#---------------------------------------------------------------------------
+
+# The QUIET tag can be used to turn on/off the messages that are generated to
+# standard output by doxygen. If QUIET is set to YES this implies that the
+# messages are off.
+# The default value is: NO.
+
+QUIET = NO
+
+# The WARNINGS tag can be used to turn on/off the warning messages that are
+# generated to standard error (stderr) by doxygen. If WARNINGS is set to YES
+# this implies that the warnings are on.
+#
+# Tip: Turn warnings on while writing the documentation.
+# The default value is: YES.
+
+WARNINGS = YES
+
+# If the WARN_IF_UNDOCUMENTED tag is set to YES then doxygen will generate
+# warnings for undocumented members. If EXTRACT_ALL is set to YES then this flag
+# will automatically be disabled.
+# The default value is: YES.
+
+WARN_IF_UNDOCUMENTED = YES
+
+# If the WARN_IF_DOC_ERROR tag is set to YES, doxygen will generate warnings for
+# potential errors in the documentation, such as not documenting some parameters
+# in a documented function, or documenting parameters that don't exist or using
+# markup commands wrongly.
+# The default value is: YES.
+
+WARN_IF_DOC_ERROR = YES
+
+# This WARN_NO_PARAMDOC option can be enabled to get warnings for functions that
+# are documented, but have no documentation for their parameters or return
+# value. If set to NO, doxygen will only warn about wrong or incomplete
+# parameter documentation, but not about the absence of documentation.
+# The default value is: NO.
+
+WARN_NO_PARAMDOC = NO
+
+# If the WARN_AS_ERROR tag is set to YES then doxygen will immediately stop when
+# a warning is encountered.
+# The default value is: NO.
+
+WARN_AS_ERROR = NO
+
+# The WARN_FORMAT tag determines the format of the warning messages that doxygen
+# can produce. The string should contain the $file, $line, and $text tags, which
+# will be replaced by the file and line number from which the warning originated
+# and the warning text. Optionally the format may contain $version, which will
+# be replaced by the version of the file (if it could be obtained via
+# FILE_VERSION_FILTER)
+# The default value is: $file:$line: $text.
+
+WARN_FORMAT = "$file:$line: $text"
+
+# The WARN_LOGFILE tag can be used to specify a file to which warning and error
+# messages should be written. If left blank the output is written to standard
+# error (stderr).
+
+WARN_LOGFILE =
+
+#---------------------------------------------------------------------------
+# Configuration options related to the input files
+#---------------------------------------------------------------------------
+
+# The INPUT tag is used to specify the files and/or directories that contain
+# documented source files. You may enter file names like myfile.cpp or
+# directories like /usr/src/myproject. Separate the files or directories with
+# spaces. See also FILE_PATTERNS and EXTENSION_MAPPING
+# Note: If this tag is empty the current directory is searched.
+
+INPUT = ../src/ ../include/RHEJ/
+
+# This tag can be used to specify the character encoding of the source files
+# that doxygen parses. Internally doxygen uses the UTF-8 encoding. Doxygen uses
+# libiconv (or the iconv built into libc) for the transcoding. See the libiconv
+# documentation (see: http://www.gnu.org/software/libiconv) for the list of
+# possible encodings.
+# The default value is: UTF-8.
+
+INPUT_ENCODING = UTF-8
+
+# If the value of the INPUT tag contains directories, you can use the
+# FILE_PATTERNS tag to specify one or more wildcard patterns (like *.cpp and
+# *.h) to filter out the source-files in the directories.
+#
+# Note that for custom extensions or not directly supported extensions you also
+# need to set EXTENSION_MAPPING for the extension otherwise the files are not
+# read by doxygen.
+#
+# If left blank the following patterns are tested:*.c, *.cc, *.cxx, *.cpp,
+# *.c++, *.java, *.ii, *.ixx, *.ipp, *.i++, *.inl, *.idl, *.ddl, *.odl, *.h,
+# *.hh, *.hxx, *.hpp, *.h++, *.cs, *.d, *.php, *.php4, *.php5, *.phtml, *.inc,
+# *.m, *.markdown, *.md, *.mm, *.dox, *.py, *.pyw, *.f90, *.f95, *.f03, *.f08,
+# *.f, *.for, *.tcl, *.vhd, *.vhdl, *.ucf and *.qsf.
+
+FILE_PATTERNS = *.c \
+ *.cc \
+ *.cxx \
+ *.cpp \
+ *.c++ \
+ *.java \
+ *.ii \
+ *.ixx \
+ *.ipp \
+ *.i++ \
+ *.inl \
+ *.idl \
+ *.ddl \
+ *.odl \
+ *.h \
+ *.hh \
+ *.hxx \
+ *.hpp \
+ *.h++ \
+ *.cs \
+ *.d \
+ *.php \
+ *.php4 \
+ *.php5 \
+ *.phtml \
+ *.inc \
+ *.m \
+ *.markdown \
+ *.md \
+ *.mm \
+ *.dox \
+ *.py \
+ *.pyw \
+ *.f90 \
+ *.f95 \
+ *.f03 \
+ *.f08 \
+ *.f \
+ *.for \
+ *.tcl \
+ *.vhd \
+ *.vhdl \
+ *.ucf \
+ *.qsf
+
+# The RECURSIVE tag can be used to specify whether or not subdirectories should
+# be searched for input files as well.
+# The default value is: NO.
+
+RECURSIVE = NO
+
+# The EXCLUDE tag can be used to specify files and/or directories that should be
+# excluded from the INPUT source files. This way you can easily exclude a
+# subdirectory from a directory tree whose root is specified with the INPUT tag.
+#
+# Note that relative paths are relative to the directory from which doxygen is
+# run.
+
+EXCLUDE =
+
+# The EXCLUDE_SYMLINKS tag can be used to select whether or not files or
+# directories that are symbolic links (a Unix file system feature) are excluded
+# from the input.
+# The default value is: NO.
+
+EXCLUDE_SYMLINKS = NO
+
+# If the value of the INPUT tag contains directories, you can use the
+# EXCLUDE_PATTERNS tag to specify one or more wildcard patterns to exclude
+# certain files from those directories.
+#
+# Note that the wildcards are matched against the file with absolute path, so to
+# exclude all test directories for example use the pattern */test/*
+
+EXCLUDE_PATTERNS =
+
+# The EXCLUDE_SYMBOLS tag can be used to specify one or more symbol names
+# (namespaces, classes, functions, etc.) that should be excluded from the
+# output. The symbol name can be a fully qualified name, a word, or if the
+# wildcard * is used, a substring. Examples: ANamespace, AClass,
+# AClass::ANamespace, ANamespace::*Test
+#
+# Note that the wildcards are matched against the file with absolute path, so to
+# exclude all test directories use the pattern */test/*
+
+EXCLUDE_SYMBOLS =
+
+# The EXAMPLE_PATH tag can be used to specify one or more files or directories
+# that contain example code fragments that are included (see the \include
+# command).
+
+EXAMPLE_PATH =
+
+# If the value of the EXAMPLE_PATH tag contains directories, you can use the
+# EXAMPLE_PATTERNS tag to specify one or more wildcard pattern (like *.cpp and
+# *.h) to filter out the source-files in the directories. If left blank all
+# files are included.
+
+EXAMPLE_PATTERNS = *
+
+# If the EXAMPLE_RECURSIVE tag is set to YES then subdirectories will be
+# searched for input files to be used with the \include or \dontinclude commands
+# irrespective of the value of the RECURSIVE tag.
+# The default value is: NO.
+
+EXAMPLE_RECURSIVE = NO
+
+# The IMAGE_PATH tag can be used to specify one or more files or directories
+# that contain images that are to be included in the documentation (see the
+# \image command).
+
+IMAGE_PATH =
+
+# The INPUT_FILTER tag can be used to specify a program that doxygen should
+# invoke to filter for each input file. Doxygen will invoke the filter program
+# by executing (via popen()) the command:
+#
+# <filter> <input-file>
+#
+# where <filter> is the value of the INPUT_FILTER tag, and <input-file> is the
+# name of an input file. Doxygen will then use the output that the filter
+# program writes to standard output. If FILTER_PATTERNS is specified, this tag
+# will be ignored.
+#
+# Note that the filter must not add or remove lines; it is applied before the
+# code is scanned, but not when the output code is generated. If lines are added
+# or removed, the anchors will not be placed correctly.
+#
+# Note that for custom extensions or not directly supported extensions you also
+# need to set EXTENSION_MAPPING for the extension otherwise the files are not
+# properly processed by doxygen.
+
+INPUT_FILTER =
+
+# The FILTER_PATTERNS tag can be used to specify filters on a per file pattern
+# basis. Doxygen will compare the file name with each pattern and apply the
+# filter if there is a match. The filters are a list of the form: pattern=filter
+# (like *.cpp=my_cpp_filter). See INPUT_FILTER for further information on how
+# filters are used. If the FILTER_PATTERNS tag is empty or if none of the
+# patterns match the file name, INPUT_FILTER is applied.
+#
+# Note that for custom extensions or not directly supported extensions you also
+# need to set EXTENSION_MAPPING for the extension otherwise the files are not
+# properly processed by doxygen.
+
+FILTER_PATTERNS =
+
+# If the FILTER_SOURCE_FILES tag is set to YES, the input filter (if set using
+# INPUT_FILTER) will also be used to filter the input files that are used for
+# producing the source files to browse (i.e. when SOURCE_BROWSER is set to YES).
+# The default value is: NO.
+
+FILTER_SOURCE_FILES = NO
+
+# The FILTER_SOURCE_PATTERNS tag can be used to specify source filters per file
+# pattern. A pattern will override the setting for FILTER_PATTERN (if any) and
+# it is also possible to disable source filtering for a specific pattern using
+# *.ext= (so without naming a filter).
+# This tag requires that the tag FILTER_SOURCE_FILES is set to YES.
+
+FILTER_SOURCE_PATTERNS =
+
+# If the USE_MDFILE_AS_MAINPAGE tag refers to the name of a markdown file that
+# is part of the input, its contents will be placed on the main page
+# (index.html). This can be useful if you have a project on for instance GitHub
+# and want to reuse the introduction page also for the doxygen output.
+
+USE_MDFILE_AS_MAINPAGE =
+
+#---------------------------------------------------------------------------
+# Configuration options related to source browsing
+#---------------------------------------------------------------------------
+
+# If the SOURCE_BROWSER tag is set to YES then a list of source files will be
+# generated. Documented entities will be cross-referenced with these sources.
+#
+# Note: To get rid of all source code in the generated output, make sure that
+# also VERBATIM_HEADERS is set to NO.
+# The default value is: NO.
+
+SOURCE_BROWSER = NO
+
+# Setting the INLINE_SOURCES tag to YES will include the body of functions,
+# classes and enums directly into the documentation.
+# The default value is: NO.
+
+INLINE_SOURCES = NO
+
+# Setting the STRIP_CODE_COMMENTS tag to YES will instruct doxygen to hide any
+# special comment blocks from generated source code fragments. Normal C, C++ and
+# Fortran comments will always remain visible.
+# The default value is: YES.
+
+STRIP_CODE_COMMENTS = YES
+
+# If the REFERENCED_BY_RELATION tag is set to YES then for each documented
+# function all documented functions referencing it will be listed.
+# The default value is: NO.
+
+REFERENCED_BY_RELATION = NO
+
+# If the REFERENCES_RELATION tag is set to YES then for each documented function
+# all documented entities called/used by that function will be listed.
+# The default value is: NO.
+
+REFERENCES_RELATION = NO
+
+# If the REFERENCES_LINK_SOURCE tag is set to YES and SOURCE_BROWSER tag is set
+# to YES then the hyperlinks from functions in REFERENCES_RELATION and
+# REFERENCED_BY_RELATION lists will link to the source code. Otherwise they will
+# link to the documentation.
+# The default value is: YES.
+
+REFERENCES_LINK_SOURCE = YES
+
+# If SOURCE_TOOLTIPS is enabled (the default) then hovering a hyperlink in the
+# source code will show a tooltip with additional information such as prototype,
+# brief description and links to the definition and documentation. Since this
+# will make the HTML file larger and loading of large files a bit slower, you
+# can opt to disable this feature.
+# The default value is: YES.
+# This tag requires that the tag SOURCE_BROWSER is set to YES.
+
+SOURCE_TOOLTIPS = YES
+
+# If the USE_HTAGS tag is set to YES then the references to source code will
+# point to the HTML generated by the htags(1) tool instead of doxygen built-in
+# source browser. The htags tool is part of GNU's global source tagging system
+# (see http://www.gnu.org/software/global/global.html). You will need version
+# 4.8.6 or higher.
+#
+# To use it do the following:
+# - Install the latest version of global
+# - Enable SOURCE_BROWSER and USE_HTAGS in the config file
+# - Make sure the INPUT points to the root of the source tree
+# - Run doxygen as normal
+#
+# Doxygen will invoke htags (and that will in turn invoke gtags), so these
+# tools must be available from the command line (i.e. in the search path).
+#
+# The result: instead of the source browser generated by doxygen, the links to
+# source code will now point to the output of htags.
+# The default value is: NO.
+# This tag requires that the tag SOURCE_BROWSER is set to YES.
+
+USE_HTAGS = NO
+
+# If the VERBATIM_HEADERS tag is set the YES then doxygen will generate a
+# verbatim copy of the header file for each class for which an include is
+# specified. Set to NO to disable this.
+# See also: Section \class.
+# The default value is: YES.
+
+VERBATIM_HEADERS = YES
+
+# If the CLANG_ASSISTED_PARSING tag is set to YES then doxygen will use the
+# clang parser (see: http://clang.llvm.org/) for more accurate parsing at the
+# cost of reduced performance. This can be particularly helpful with template
+# rich C++ code for which doxygen's built-in parser lacks the necessary type
+# information.
+# Note: The availability of this option depends on whether or not doxygen was
+# generated with the -Duse-libclang=ON option for CMake.
+# The default value is: NO.
+
+CLANG_ASSISTED_PARSING = NO
+
+# If clang assisted parsing is enabled you can provide the compiler with command
+# line options that you would normally use when invoking the compiler. Note that
+# the include paths will already be set by doxygen for the files and directories
+# specified with INPUT and INCLUDE_PATH.
+# This tag requires that the tag CLANG_ASSISTED_PARSING is set to YES.
+
+CLANG_OPTIONS =
+
+#---------------------------------------------------------------------------
+# Configuration options related to the alphabetical class index
+#---------------------------------------------------------------------------
+
+# If the ALPHABETICAL_INDEX tag is set to YES, an alphabetical index of all
+# compounds will be generated. Enable this if the project contains a lot of
+# classes, structs, unions or interfaces.
+# The default value is: YES.
+
+ALPHABETICAL_INDEX = YES
+
+# The COLS_IN_ALPHA_INDEX tag can be used to specify the number of columns in
+# which the alphabetical index list will be split.
+# Minimum value: 1, maximum value: 20, default value: 5.
+# This tag requires that the tag ALPHABETICAL_INDEX is set to YES.
+
+COLS_IN_ALPHA_INDEX = 5
+
+# In case all classes in a project start with a common prefix, all classes will
+# be put under the same header in the alphabetical index. The IGNORE_PREFIX tag
+# can be used to specify a prefix (or a list of prefixes) that should be ignored
+# while generating the index headers.
+# This tag requires that the tag ALPHABETICAL_INDEX is set to YES.
+
+IGNORE_PREFIX =
+
+#---------------------------------------------------------------------------
+# Configuration options related to the HTML output
+#---------------------------------------------------------------------------
+
+# If the GENERATE_HTML tag is set to YES, doxygen will generate HTML output
+# The default value is: YES.
+
+GENERATE_HTML = YES
+
+# The HTML_OUTPUT tag is used to specify where the HTML docs will be put. If a
+# relative path is entered the value of OUTPUT_DIRECTORY will be put in front of
+# it.
+# The default directory is: html.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+HTML_OUTPUT = html
+
+# The HTML_FILE_EXTENSION tag can be used to specify the file extension for each
+# generated HTML page (for example: .htm, .php, .asp).
+# The default value is: .html.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+HTML_FILE_EXTENSION = .html
+
+# The HTML_HEADER tag can be used to specify a user-defined HTML header file for
+# each generated HTML page. If the tag is left blank doxygen will generate a
+# standard header.
+#
+# To get valid HTML the header file that includes any scripts and style sheets
+# that doxygen needs, which is dependent on the configuration options used (e.g.
+# the setting GENERATE_TREEVIEW). It is highly recommended to start with a
+# default header using
+# doxygen -w html new_header.html new_footer.html new_stylesheet.css
+# YourConfigFile
+# and then modify the file new_header.html. See also section "Doxygen usage"
+# for information on how to generate the default header that doxygen normally
+# uses.
+# Note: The header is subject to change so you typically have to regenerate the
+# default header when upgrading to a newer version of doxygen. For a description
+# of the possible markers and block names see the documentation.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+HTML_HEADER =
+
+# The HTML_FOOTER tag can be used to specify a user-defined HTML footer for each
+# generated HTML page. If the tag is left blank doxygen will generate a standard
+# footer. See HTML_HEADER for more information on how to generate a default
+# footer and what special commands can be used inside the footer. See also
+# section "Doxygen usage" for information on how to generate the default footer
+# that doxygen normally uses.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+HTML_FOOTER =
+
+# The HTML_STYLESHEET tag can be used to specify a user-defined cascading style
+# sheet that is used by each HTML page. It can be used to fine-tune the look of
+# the HTML output. If left blank doxygen will generate a default style sheet.
+# See also section "Doxygen usage" for information on how to generate the style
+# sheet that doxygen normally uses.
+# Note: It is recommended to use HTML_EXTRA_STYLESHEET instead of this tag, as
+# it is more robust and this tag (HTML_STYLESHEET) will in the future become
+# obsolete.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+HTML_STYLESHEET =
+
+# The HTML_EXTRA_STYLESHEET tag can be used to specify additional user-defined
+# cascading style sheets that are included after the standard style sheets
+# created by doxygen. Using this option one can overrule certain style aspects.
+# This is preferred over using HTML_STYLESHEET since it does not replace the
+# standard style sheet and is therefore more robust against future updates.
+# Doxygen will copy the style sheet files to the output directory.
+# Note: The order of the extra style sheet files is of importance (e.g. the last
+# style sheet in the list overrules the setting of the previous ones in the
+# list). For an example see the documentation.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+HTML_EXTRA_STYLESHEET =
+
+# The HTML_EXTRA_FILES tag can be used to specify one or more extra images or
+# other source files which should be copied to the HTML output directory. Note
+# that these files will be copied to the base HTML output directory. Use the
+# $relpath^ marker in the HTML_HEADER and/or HTML_FOOTER files to load these
+# files. In the HTML_STYLESHEET file, use the file name only. Also note that the
+# files will be copied as-is; there are no commands or markers available.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+HTML_EXTRA_FILES =
+
+# The HTML_COLORSTYLE_HUE tag controls the color of the HTML output. Doxygen
+# will adjust the colors in the style sheet and background images according to
+# this color. Hue is specified as an angle on a colorwheel, see
+# http://en.wikipedia.org/wiki/Hue for more information. For instance the value
+# 0 represents red, 60 is yellow, 120 is green, 180 is cyan, 240 is blue, 300
+# purple, and 360 is red again.
+# Minimum value: 0, maximum value: 359, default value: 220.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+HTML_COLORSTYLE_HUE = 220
+
+# The HTML_COLORSTYLE_SAT tag controls the purity (or saturation) of the colors
+# in the HTML output. For a value of 0 the output will use grayscales only. A
+# value of 255 will produce the most vivid colors.
+# Minimum value: 0, maximum value: 255, default value: 100.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+HTML_COLORSTYLE_SAT = 100
+
+# The HTML_COLORSTYLE_GAMMA tag controls the gamma correction applied to the
+# luminance component of the colors in the HTML output. Values below 100
+# gradually make the output lighter, whereas values above 100 make the output
+# darker. The value divided by 100 is the actual gamma applied, so 80 represents
+# a gamma of 0.8, The value 220 represents a gamma of 2.2, and 100 does not
+# change the gamma.
+# Minimum value: 40, maximum value: 240, default value: 80.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+HTML_COLORSTYLE_GAMMA = 80
+
+# If the HTML_TIMESTAMP tag is set to YES then the footer of each generated HTML
+# page will contain the date and time when the page was generated. Setting this
+# to YES can help to show when doxygen was last run and thus if the
+# documentation is up to date.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+HTML_TIMESTAMP = NO
+
+# If the HTML_DYNAMIC_SECTIONS tag is set to YES then the generated HTML
+# documentation will contain sections that can be hidden and shown after the
+# page has loaded.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+HTML_DYNAMIC_SECTIONS = NO
+
+# With HTML_INDEX_NUM_ENTRIES one can control the preferred number of entries
+# shown in the various tree structured indices initially; the user can expand
+# and collapse entries dynamically later on. Doxygen will expand the tree to
+# such a level that at most the specified number of entries are visible (unless
+# a fully collapsed tree already exceeds this amount). So setting the number of
+# entries 1 will produce a full collapsed tree by default. 0 is a special value
+# representing an infinite number of entries and will result in a full expanded
+# tree by default.
+# Minimum value: 0, maximum value: 9999, default value: 100.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+HTML_INDEX_NUM_ENTRIES = 100
+
+# If the GENERATE_DOCSET tag is set to YES, additional index files will be
+# generated that can be used as input for Apple's Xcode 3 integrated development
+# environment (see: http://developer.apple.com/tools/xcode/), introduced with
+# OSX 10.5 (Leopard). To create a documentation set, doxygen will generate a
+# Makefile in the HTML output directory. Running make will produce the docset in
+# that directory and running make install will install the docset in
+# ~/Library/Developer/Shared/Documentation/DocSets so that Xcode will find it at
+# startup. See http://developer.apple.com/tools/creatingdocsetswithdoxygen.html
+# for more information.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+GENERATE_DOCSET = NO
+
+# This tag determines the name of the docset feed. A documentation feed provides
+# an umbrella under which multiple documentation sets from a single provider
+# (such as a company or product suite) can be grouped.
+# The default value is: Doxygen generated docs.
+# This tag requires that the tag GENERATE_DOCSET is set to YES.
+
+DOCSET_FEEDNAME = "Doxygen generated docs"
+
+# This tag specifies a string that should uniquely identify the documentation
+# set bundle. This should be a reverse domain-name style string, e.g.
+# com.mycompany.MyDocSet. Doxygen will append .docset to the name.
+# The default value is: org.doxygen.Project.
+# This tag requires that the tag GENERATE_DOCSET is set to YES.
+
+DOCSET_BUNDLE_ID = org.doxygen.Project
+
+# The DOCSET_PUBLISHER_ID tag specifies a string that should uniquely identify
+# the documentation publisher. This should be a reverse domain-name style
+# string, e.g. com.mycompany.MyDocSet.documentation.
+# The default value is: org.doxygen.Publisher.
+# This tag requires that the tag GENERATE_DOCSET is set to YES.
+
+DOCSET_PUBLISHER_ID = org.doxygen.Publisher
+
+# The DOCSET_PUBLISHER_NAME tag identifies the documentation publisher.
+# The default value is: Publisher.
+# This tag requires that the tag GENERATE_DOCSET is set to YES.
+
+DOCSET_PUBLISHER_NAME = Publisher
+
+# If the GENERATE_HTMLHELP tag is set to YES then doxygen generates three
+# additional HTML index files: index.hhp, index.hhc, and index.hhk. The
+# index.hhp is a project file that can be read by Microsoft's HTML Help Workshop
+# (see: http://www.microsoft.com/en-us/download/details.aspx?id=21138) on
+# Windows.
+#
+# The HTML Help Workshop contains a compiler that can convert all HTML output
+# generated by doxygen into a single compiled HTML file (.chm). Compiled HTML
+# files are now used as the Windows 98 help format, and will replace the old
+# Windows help format (.hlp) on all Windows platforms in the future. Compressed
+# HTML files also contain an index, a table of contents, and you can search for
+# words in the documentation. The HTML workshop also contains a viewer for
+# compressed HTML files.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+GENERATE_HTMLHELP = NO
+
+# The CHM_FILE tag can be used to specify the file name of the resulting .chm
+# file. You can add a path in front of the file if the result should not be
+# written to the html output directory.
+# This tag requires that the tag GENERATE_HTMLHELP is set to YES.
+
+CHM_FILE =
+
+# The HHC_LOCATION tag can be used to specify the location (absolute path
+# including file name) of the HTML help compiler (hhc.exe). If non-empty,
+# doxygen will try to run the HTML help compiler on the generated index.hhp.
+# The file has to be specified with full path.
+# This tag requires that the tag GENERATE_HTMLHELP is set to YES.
+
+HHC_LOCATION =
+
+# The GENERATE_CHI flag controls if a separate .chi index file is generated
+# (YES) or that it should be included in the master .chm file (NO).
+# The default value is: NO.
+# This tag requires that the tag GENERATE_HTMLHELP is set to YES.
+
+GENERATE_CHI = NO
+
+# The CHM_INDEX_ENCODING is used to encode HtmlHelp index (hhk), content (hhc)
+# and project file content.
+# This tag requires that the tag GENERATE_HTMLHELP is set to YES.
+
+CHM_INDEX_ENCODING =
+
+# The BINARY_TOC flag controls whether a binary table of contents is generated
+# (YES) or a normal table of contents (NO) in the .chm file. Furthermore it
+# enables the Previous and Next buttons.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_HTMLHELP is set to YES.
+
+BINARY_TOC = NO
+
+# The TOC_EXPAND flag can be set to YES to add extra items for group members to
+# the table of contents of the HTML help documentation and to the tree view.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_HTMLHELP is set to YES.
+
+TOC_EXPAND = NO
+
+# If the GENERATE_QHP tag is set to YES and both QHP_NAMESPACE and
+# QHP_VIRTUAL_FOLDER are set, an additional index file will be generated that
+# can be used as input for Qt's qhelpgenerator to generate a Qt Compressed Help
+# (.qch) of the generated HTML documentation.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+GENERATE_QHP = NO
+
+# If the QHG_LOCATION tag is specified, the QCH_FILE tag can be used to specify
+# the file name of the resulting .qch file. The path specified is relative to
+# the HTML output folder.
+# This tag requires that the tag GENERATE_QHP is set to YES.
+
+QCH_FILE =
+
+# The QHP_NAMESPACE tag specifies the namespace to use when generating Qt Help
+# Project output. For more information please see Qt Help Project / Namespace
+# (see: http://qt-project.org/doc/qt-4.8/qthelpproject.html#namespace).
+# The default value is: org.doxygen.Project.
+# This tag requires that the tag GENERATE_QHP is set to YES.
+
+QHP_NAMESPACE = org.doxygen.Project
+
+# The QHP_VIRTUAL_FOLDER tag specifies the namespace to use when generating Qt
+# Help Project output. For more information please see Qt Help Project / Virtual
+# Folders (see: http://qt-project.org/doc/qt-4.8/qthelpproject.html#virtual-
+# folders).
+# The default value is: doc.
+# This tag requires that the tag GENERATE_QHP is set to YES.
+
+QHP_VIRTUAL_FOLDER = doc
+
+# If the QHP_CUST_FILTER_NAME tag is set, it specifies the name of a custom
+# filter to add. For more information please see Qt Help Project / Custom
+# Filters (see: http://qt-project.org/doc/qt-4.8/qthelpproject.html#custom-
+# filters).
+# This tag requires that the tag GENERATE_QHP is set to YES.
+
+QHP_CUST_FILTER_NAME =
+
+# The QHP_CUST_FILTER_ATTRS tag specifies the list of the attributes of the
+# custom filter to add. For more information please see Qt Help Project / Custom
+# Filters (see: http://qt-project.org/doc/qt-4.8/qthelpproject.html#custom-
+# filters).
+# This tag requires that the tag GENERATE_QHP is set to YES.
+
+QHP_CUST_FILTER_ATTRS =
+
+# The QHP_SECT_FILTER_ATTRS tag specifies the list of the attributes this
+# project's filter section matches. Qt Help Project / Filter Attributes (see:
+# http://qt-project.org/doc/qt-4.8/qthelpproject.html#filter-attributes).
+# This tag requires that the tag GENERATE_QHP is set to YES.
+
+QHP_SECT_FILTER_ATTRS =
+
+# The QHG_LOCATION tag can be used to specify the location of Qt's
+# qhelpgenerator. If non-empty doxygen will try to run qhelpgenerator on the
+# generated .qhp file.
+# This tag requires that the tag GENERATE_QHP is set to YES.
+
+QHG_LOCATION =
+
+# If the GENERATE_ECLIPSEHELP tag is set to YES, additional index files will be
+# generated, together with the HTML files, they form an Eclipse help plugin. To
+# install this plugin and make it available under the help contents menu in
+# Eclipse, the contents of the directory containing the HTML and XML files needs
+# to be copied into the plugins directory of eclipse. The name of the directory
+# within the plugins directory should be the same as the ECLIPSE_DOC_ID value.
+# After copying Eclipse needs to be restarted before the help appears.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+GENERATE_ECLIPSEHELP = NO
+
+# A unique identifier for the Eclipse help plugin. When installing the plugin
+# the directory name containing the HTML and XML files should also have this
+# name. Each documentation set should have its own identifier.
+# The default value is: org.doxygen.Project.
+# This tag requires that the tag GENERATE_ECLIPSEHELP is set to YES.
+
+ECLIPSE_DOC_ID = org.doxygen.Project
+
+# If you want full control over the layout of the generated HTML pages it might
+# be necessary to disable the index and replace it with your own. The
+# DISABLE_INDEX tag can be used to turn on/off the condensed index (tabs) at top
+# of each HTML page. A value of NO enables the index and the value YES disables
+# it. Since the tabs in the index contain the same information as the navigation
+# tree, you can set this option to YES if you also set GENERATE_TREEVIEW to YES.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+DISABLE_INDEX = NO
+
+# The GENERATE_TREEVIEW tag is used to specify whether a tree-like index
+# structure should be generated to display hierarchical information. If the tag
+# value is set to YES, a side panel will be generated containing a tree-like
+# index structure (just like the one that is generated for HTML Help). For this
+# to work a browser that supports JavaScript, DHTML, CSS and frames is required
+# (i.e. any modern browser). Windows users are probably better off using the
+# HTML help feature. Via custom style sheets (see HTML_EXTRA_STYLESHEET) one can
+# further fine-tune the look of the index. As an example, the default style
+# sheet generated by doxygen has an example that shows how to put an image at
+# the root of the tree instead of the PROJECT_NAME. Since the tree basically has
+# the same information as the tab index, you could consider setting
+# DISABLE_INDEX to YES when enabling this option.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+GENERATE_TREEVIEW = NO
+
+# The ENUM_VALUES_PER_LINE tag can be used to set the number of enum values that
+# doxygen will group on one line in the generated HTML documentation.
+#
+# Note that a value of 0 will completely suppress the enum values from appearing
+# in the overview section.
+# Minimum value: 0, maximum value: 20, default value: 4.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+ENUM_VALUES_PER_LINE = 4
+
+# If the treeview is enabled (see GENERATE_TREEVIEW) then this tag can be used
+# to set the initial width (in pixels) of the frame in which the tree is shown.
+# Minimum value: 0, maximum value: 1500, default value: 250.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+TREEVIEW_WIDTH = 250
+
+# If the EXT_LINKS_IN_WINDOW option is set to YES, doxygen will open links to
+# external symbols imported via tag files in a separate window.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+EXT_LINKS_IN_WINDOW = NO
+
+# Use this tag to change the font size of LaTeX formulas included as images in
+# the HTML documentation. When you change the font size after a successful
+# doxygen run you need to manually remove any form_*.png images from the HTML
+# output directory to force them to be regenerated.
+# Minimum value: 8, maximum value: 50, default value: 10.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+FORMULA_FONTSIZE = 10
+
+# Use the FORMULA_TRANPARENT tag to determine whether or not the images
+# generated for formulas are transparent PNGs. Transparent PNGs are not
+# supported properly for IE 6.0, but are supported on all modern browsers.
+#
+# Note that when changing this option you need to delete any form_*.png files in
+# the HTML output directory before the changes have effect.
+# The default value is: YES.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+FORMULA_TRANSPARENT = YES
+
+# Enable the USE_MATHJAX option to render LaTeX formulas using MathJax (see
+# http://www.mathjax.org) which uses client side Javascript for the rendering
+# instead of using pre-rendered bitmaps. Use this if you do not have LaTeX
+# installed or if you want to formulas look prettier in the HTML output. When
+# enabled you may also need to install MathJax separately and configure the path
+# to it using the MATHJAX_RELPATH option.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+USE_MATHJAX = NO
+
+# When MathJax is enabled you can set the default output format to be used for
+# the MathJax output. See the MathJax site (see:
+# http://docs.mathjax.org/en/latest/output.html) for more details.
+# Possible values are: HTML-CSS (which is slower, but has the best
+# compatibility), NativeMML (i.e. MathML) and SVG.
+# The default value is: HTML-CSS.
+# This tag requires that the tag USE_MATHJAX is set to YES.
+
+MATHJAX_FORMAT = HTML-CSS
+
+# When MathJax is enabled you need to specify the location relative to the HTML
+# output directory using the MATHJAX_RELPATH option. The destination directory
+# should contain the MathJax.js script. For instance, if the mathjax directory
+# is located at the same level as the HTML output directory, then
+# MATHJAX_RELPATH should be ../mathjax. The default value points to the MathJax
+# Content Delivery Network so you can quickly see the result without installing
+# MathJax. However, it is strongly recommended to install a local copy of
+# MathJax from http://www.mathjax.org before deployment.
+# The default value is: http://cdn.mathjax.org/mathjax/latest.
+# This tag requires that the tag USE_MATHJAX is set to YES.
+
+MATHJAX_RELPATH = http://cdn.mathjax.org/mathjax/latest
+
+# The MATHJAX_EXTENSIONS tag can be used to specify one or more MathJax
+# extension names that should be enabled during MathJax rendering. For example
+# MATHJAX_EXTENSIONS = TeX/AMSmath TeX/AMSsymbols
+# This tag requires that the tag USE_MATHJAX is set to YES.
+
+MATHJAX_EXTENSIONS =
+
+# The MATHJAX_CODEFILE tag can be used to specify a file with javascript pieces
+# of code that will be used on startup of the MathJax code. See the MathJax site
+# (see: http://docs.mathjax.org/en/latest/output.html) for more details. For an
+# example see the documentation.
+# This tag requires that the tag USE_MATHJAX is set to YES.
+
+MATHJAX_CODEFILE =
+
+# When the SEARCHENGINE tag is enabled doxygen will generate a search box for
+# the HTML output. The underlying search engine uses javascript and DHTML and
+# should work on any modern browser. Note that when using HTML help
+# (GENERATE_HTMLHELP), Qt help (GENERATE_QHP), or docsets (GENERATE_DOCSET)
+# there is already a search function so this one should typically be disabled.
+# For large projects the javascript based search engine can be slow, then
+# enabling SERVER_BASED_SEARCH may provide a better solution. It is possible to
+# search using the keyboard; to jump to the search box use <access key> + S
+# (what the <access key> is depends on the OS and browser, but it is typically
+# <CTRL>, <ALT>/<option>, or both). Inside the search box use the <cursor down
+# key> to jump into the search results window, the results can be navigated
+# using the <cursor keys>. Press <Enter> to select an item or <escape> to cancel
+# the search. The filter options can be selected when the cursor is inside the
+# search box by pressing <Shift>+<cursor down>. Also here use the <cursor keys>
+# to select a filter and <Enter> or <escape> to activate or cancel the filter
+# option.
+# The default value is: YES.
+# This tag requires that the tag GENERATE_HTML is set to YES.
+
+SEARCHENGINE = YES
+
+# When the SERVER_BASED_SEARCH tag is enabled the search engine will be
+# implemented using a web server instead of a web client using Javascript. There
+# are two flavors of web server based searching depending on the EXTERNAL_SEARCH
+# setting. When disabled, doxygen will generate a PHP script for searching and
+# an index file used by the script. When EXTERNAL_SEARCH is enabled the indexing
+# and searching needs to be provided by external tools. See the section
+# "External Indexing and Searching" for details.
+# The default value is: NO.
+# This tag requires that the tag SEARCHENGINE is set to YES.
+
+SERVER_BASED_SEARCH = NO
+
+# When EXTERNAL_SEARCH tag is enabled doxygen will no longer generate the PHP
+# script for searching. Instead the search results are written to an XML file
+# which needs to be processed by an external indexer. Doxygen will invoke an
+# external search engine pointed to by the SEARCHENGINE_URL option to obtain the
+# search results.
+#
+# Doxygen ships with an example indexer (doxyindexer) and search engine
+# (doxysearch.cgi) which are based on the open source search engine library
+# Xapian (see: http://xapian.org/).
+#
+# See the section "External Indexing and Searching" for details.
+# The default value is: NO.
+# This tag requires that the tag SEARCHENGINE is set to YES.
+
+EXTERNAL_SEARCH = NO
+
+# The SEARCHENGINE_URL should point to a search engine hosted by a web server
+# which will return the search results when EXTERNAL_SEARCH is enabled.
+#
+# Doxygen ships with an example indexer (doxyindexer) and search engine
+# (doxysearch.cgi) which are based on the open source search engine library
+# Xapian (see: http://xapian.org/). See the section "External Indexing and
+# Searching" for details.
+# This tag requires that the tag SEARCHENGINE is set to YES.
+
+SEARCHENGINE_URL =
+
+# When SERVER_BASED_SEARCH and EXTERNAL_SEARCH are both enabled the unindexed
+# search data is written to a file for indexing by an external tool. With the
+# SEARCHDATA_FILE tag the name of this file can be specified.
+# The default file is: searchdata.xml.
+# This tag requires that the tag SEARCHENGINE is set to YES.
+
+SEARCHDATA_FILE = searchdata.xml
+
+# When SERVER_BASED_SEARCH and EXTERNAL_SEARCH are both enabled the
+# EXTERNAL_SEARCH_ID tag can be used as an identifier for the project. This is
+# useful in combination with EXTRA_SEARCH_MAPPINGS to search through multiple
+# projects and redirect the results back to the right project.
+# This tag requires that the tag SEARCHENGINE is set to YES.
+
+EXTERNAL_SEARCH_ID =
+
+# The EXTRA_SEARCH_MAPPINGS tag can be used to enable searching through doxygen
+# projects other than the one defined by this configuration file, but that are
+# all added to the same external search index. Each project needs to have a
+# unique id set via EXTERNAL_SEARCH_ID. The search mapping then maps the id of
+# to a relative location where the documentation can be found. The format is:
+# EXTRA_SEARCH_MAPPINGS = tagname1=loc1 tagname2=loc2 ...
+# This tag requires that the tag SEARCHENGINE is set to YES.
+
+EXTRA_SEARCH_MAPPINGS =
+
+#---------------------------------------------------------------------------
+# Configuration options related to the LaTeX output
+#---------------------------------------------------------------------------
+
+# If the GENERATE_LATEX tag is set to YES, doxygen will generate LaTeX output.
+# The default value is: YES.
+
+GENERATE_LATEX = YES
+
+# The LATEX_OUTPUT tag is used to specify where the LaTeX docs will be put. If a
+# relative path is entered the value of OUTPUT_DIRECTORY will be put in front of
+# it.
+# The default directory is: latex.
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+LATEX_OUTPUT = latex
+
+# The LATEX_CMD_NAME tag can be used to specify the LaTeX command name to be
+# invoked.
+#
+# Note that when enabling USE_PDFLATEX this option is only used for generating
+# bitmaps for formulas in the HTML output, but not in the Makefile that is
+# written to the output directory.
+# The default file is: latex.
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+LATEX_CMD_NAME = latex
+
+# The MAKEINDEX_CMD_NAME tag can be used to specify the command name to generate
+# index for LaTeX.
+# The default file is: makeindex.
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+MAKEINDEX_CMD_NAME = makeindex
+
+# If the COMPACT_LATEX tag is set to YES, doxygen generates more compact LaTeX
+# documents. This may be useful for small projects and may help to save some
+# trees in general.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+COMPACT_LATEX = NO
+
+# The PAPER_TYPE tag can be used to set the paper type that is used by the
+# printer.
+# Possible values are: a4 (210 x 297 mm), letter (8.5 x 11 inches), legal (8.5 x
+# 14 inches) and executive (7.25 x 10.5 inches).
+# The default value is: a4.
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+PAPER_TYPE = a4
+
+# The EXTRA_PACKAGES tag can be used to specify one or more LaTeX package names
+# that should be included in the LaTeX output. The package can be specified just
+# by its name or with the correct syntax as to be used with the LaTeX
+# \usepackage command. To get the times font for instance you can specify :
+# EXTRA_PACKAGES=times or EXTRA_PACKAGES={times}
+# To use the option intlimits with the amsmath package you can specify:
+# EXTRA_PACKAGES=[intlimits]{amsmath}
+# If left blank no extra packages will be included.
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+EXTRA_PACKAGES =
+
+# The LATEX_HEADER tag can be used to specify a personal LaTeX header for the
+# generated LaTeX document. The header should contain everything until the first
+# chapter. If it is left blank doxygen will generate a standard header. See
+# section "Doxygen usage" for information on how to let doxygen write the
+# default header to a separate file.
+#
+# Note: Only use a user-defined header if you know what you are doing! The
+# following commands have a special meaning inside the header: $title,
+# $datetime, $date, $doxygenversion, $projectname, $projectnumber,
+# $projectbrief, $projectlogo. Doxygen will replace $title with the empty
+# string, for the replacement values of the other commands the user is referred
+# to HTML_HEADER.
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+LATEX_HEADER =
+
+# The LATEX_FOOTER tag can be used to specify a personal LaTeX footer for the
+# generated LaTeX document. The footer should contain everything after the last
+# chapter. If it is left blank doxygen will generate a standard footer. See
+# LATEX_HEADER for more information on how to generate a default footer and what
+# special commands can be used inside the footer.
+#
+# Note: Only use a user-defined footer if you know what you are doing!
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+LATEX_FOOTER =
+
+# The LATEX_EXTRA_STYLESHEET tag can be used to specify additional user-defined
+# LaTeX style sheets that are included after the standard style sheets created
+# by doxygen. Using this option one can overrule certain style aspects. Doxygen
+# will copy the style sheet files to the output directory.
+# Note: The order of the extra style sheet files is of importance (e.g. the last
+# style sheet in the list overrules the setting of the previous ones in the
+# list).
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+LATEX_EXTRA_STYLESHEET =
+
+# The LATEX_EXTRA_FILES tag can be used to specify one or more extra images or
+# other source files which should be copied to the LATEX_OUTPUT output
+# directory. Note that the files will be copied as-is; there are no commands or
+# markers available.
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+LATEX_EXTRA_FILES =
+
+# If the PDF_HYPERLINKS tag is set to YES, the LaTeX that is generated is
+# prepared for conversion to PDF (using ps2pdf or pdflatex). The PDF file will
+# contain links (just like the HTML output) instead of page references. This
+# makes the output suitable for online browsing using a PDF viewer.
+# The default value is: YES.
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+PDF_HYPERLINKS = YES
+
+# If the USE_PDFLATEX tag is set to YES, doxygen will use pdflatex to generate
+# the PDF file directly from the LaTeX files. Set this option to YES, to get a
+# higher quality PDF documentation.
+# The default value is: YES.
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+USE_PDFLATEX = YES
+
+# If the LATEX_BATCHMODE tag is set to YES, doxygen will add the \batchmode
+# command to the generated LaTeX files. This will instruct LaTeX to keep running
+# if errors occur, instead of asking the user for help. This option is also used
+# when generating formulas in HTML.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+LATEX_BATCHMODE = NO
+
+# If the LATEX_HIDE_INDICES tag is set to YES then doxygen will not include the
+# index chapters (such as File Index, Compound Index, etc.) in the output.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+LATEX_HIDE_INDICES = NO
+
+# If the LATEX_SOURCE_CODE tag is set to YES then doxygen will include source
+# code with syntax highlighting in the LaTeX output.
+#
+# Note that which sources are shown also depends on other settings such as
+# SOURCE_BROWSER.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+LATEX_SOURCE_CODE = NO
+
+# The LATEX_BIB_STYLE tag can be used to specify the style to use for the
+# bibliography, e.g. plainnat, or ieeetr. See
+# http://en.wikipedia.org/wiki/BibTeX and \cite for more info.
+# The default value is: plain.
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+LATEX_BIB_STYLE = plain
+
+# If the LATEX_TIMESTAMP tag is set to YES then the footer of each generated
+# page will contain the date and time when the page was generated. Setting this
+# to NO can help when comparing the output of multiple runs.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_LATEX is set to YES.
+
+LATEX_TIMESTAMP = NO
+
+#---------------------------------------------------------------------------
+# Configuration options related to the RTF output
+#---------------------------------------------------------------------------
+
+# If the GENERATE_RTF tag is set to YES, doxygen will generate RTF output. The
+# RTF output is optimized for Word 97 and may not look too pretty with other RTF
+# readers/editors.
+# The default value is: NO.
+
+GENERATE_RTF = NO
+
+# The RTF_OUTPUT tag is used to specify where the RTF docs will be put. If a
+# relative path is entered the value of OUTPUT_DIRECTORY will be put in front of
+# it.
+# The default directory is: rtf.
+# This tag requires that the tag GENERATE_RTF is set to YES.
+
+RTF_OUTPUT = rtf
+
+# If the COMPACT_RTF tag is set to YES, doxygen generates more compact RTF
+# documents. This may be useful for small projects and may help to save some
+# trees in general.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_RTF is set to YES.
+
+COMPACT_RTF = NO
+
+# If the RTF_HYPERLINKS tag is set to YES, the RTF that is generated will
+# contain hyperlink fields. The RTF file will contain links (just like the HTML
+# output) instead of page references. This makes the output suitable for online
+# browsing using Word or some other Word compatible readers that support those
+# fields.
+#
+# Note: WordPad (write) and others do not support links.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_RTF is set to YES.
+
+RTF_HYPERLINKS = NO
+
+# Load stylesheet definitions from file. Syntax is similar to doxygen's config
+# file, i.e. a series of assignments. You only have to provide replacements,
+# missing definitions are set to their default value.
+#
+# See also section "Doxygen usage" for information on how to generate the
+# default style sheet that doxygen normally uses.
+# This tag requires that the tag GENERATE_RTF is set to YES.
+
+RTF_STYLESHEET_FILE =
+
+# Set optional variables used in the generation of an RTF document. Syntax is
+# similar to doxygen's config file. A template extensions file can be generated
+# using doxygen -e rtf extensionFile.
+# This tag requires that the tag GENERATE_RTF is set to YES.
+
+RTF_EXTENSIONS_FILE =
+
+# If the RTF_SOURCE_CODE tag is set to YES then doxygen will include source code
+# with syntax highlighting in the RTF output.
+#
+# Note that which sources are shown also depends on other settings such as
+# SOURCE_BROWSER.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_RTF is set to YES.
+
+RTF_SOURCE_CODE = NO
+
+#---------------------------------------------------------------------------
+# Configuration options related to the man page output
+#---------------------------------------------------------------------------
+
+# If the GENERATE_MAN tag is set to YES, doxygen will generate man pages for
+# classes and files.
+# The default value is: NO.
+
+GENERATE_MAN = NO
+
+# The MAN_OUTPUT tag is used to specify where the man pages will be put. If a
+# relative path is entered the value of OUTPUT_DIRECTORY will be put in front of
+# it. A directory man3 will be created inside the directory specified by
+# MAN_OUTPUT.
+# The default directory is: man.
+# This tag requires that the tag GENERATE_MAN is set to YES.
+
+MAN_OUTPUT = man
+
+# The MAN_EXTENSION tag determines the extension that is added to the generated
+# man pages. In case the manual section does not start with a number, the number
+# 3 is prepended. The dot (.) at the beginning of the MAN_EXTENSION tag is
+# optional.
+# The default value is: .3.
+# This tag requires that the tag GENERATE_MAN is set to YES.
+
+MAN_EXTENSION = .3
+
+# The MAN_SUBDIR tag determines the name of the directory created within
+# MAN_OUTPUT in which the man pages are placed. If defaults to man followed by
+# MAN_EXTENSION with the initial . removed.
+# This tag requires that the tag GENERATE_MAN is set to YES.
+
+MAN_SUBDIR =
+
+# If the MAN_LINKS tag is set to YES and doxygen generates man output, then it
+# will generate one additional man file for each entity documented in the real
+# man page(s). These additional files only source the real man page, but without
+# them the man command would be unable to find the correct page.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_MAN is set to YES.
+
+MAN_LINKS = NO
+
+#---------------------------------------------------------------------------
+# Configuration options related to the XML output
+#---------------------------------------------------------------------------
+
+# If the GENERATE_XML tag is set to YES, doxygen will generate an XML file that
+# captures the structure of the code including all documentation.
+# The default value is: NO.
+
+GENERATE_XML = NO
+
+# The XML_OUTPUT tag is used to specify where the XML pages will be put. If a
+# relative path is entered the value of OUTPUT_DIRECTORY will be put in front of
+# it.
+# The default directory is: xml.
+# This tag requires that the tag GENERATE_XML is set to YES.
+
+XML_OUTPUT = xml
+
+# If the XML_PROGRAMLISTING tag is set to YES, doxygen will dump the program
+# listings (including syntax highlighting and cross-referencing information) to
+# the XML output. Note that enabling this will significantly increase the size
+# of the XML output.
+# The default value is: YES.
+# This tag requires that the tag GENERATE_XML is set to YES.
+
+XML_PROGRAMLISTING = YES
+
+#---------------------------------------------------------------------------
+# Configuration options related to the DOCBOOK output
+#---------------------------------------------------------------------------
+
+# If the GENERATE_DOCBOOK tag is set to YES, doxygen will generate Docbook files
+# that can be used to generate PDF.
+# The default value is: NO.
+
+GENERATE_DOCBOOK = NO
+
+# The DOCBOOK_OUTPUT tag is used to specify where the Docbook pages will be put.
+# If a relative path is entered the value of OUTPUT_DIRECTORY will be put in
+# front of it.
+# The default directory is: docbook.
+# This tag requires that the tag GENERATE_DOCBOOK is set to YES.
+
+DOCBOOK_OUTPUT = docbook
+
+# If the DOCBOOK_PROGRAMLISTING tag is set to YES, doxygen will include the
+# program listings (including syntax highlighting and cross-referencing
+# information) to the DOCBOOK output. Note that enabling this will significantly
+# increase the size of the DOCBOOK output.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_DOCBOOK is set to YES.
+
+DOCBOOK_PROGRAMLISTING = NO
+
+#---------------------------------------------------------------------------
+# Configuration options for the AutoGen Definitions output
+#---------------------------------------------------------------------------
+
+# If the GENERATE_AUTOGEN_DEF tag is set to YES, doxygen will generate an
+# AutoGen Definitions (see http://autogen.sf.net) file that captures the
+# structure of the code including all documentation. Note that this feature is
+# still experimental and incomplete at the moment.
+# The default value is: NO.
+
+GENERATE_AUTOGEN_DEF = NO
+
+#---------------------------------------------------------------------------
+# Configuration options related to the Perl module output
+#---------------------------------------------------------------------------
+
+# If the GENERATE_PERLMOD tag is set to YES, doxygen will generate a Perl module
+# file that captures the structure of the code including all documentation.
+#
+# Note that this feature is still experimental and incomplete at the moment.
+# The default value is: NO.
+
+GENERATE_PERLMOD = NO
+
+# If the PERLMOD_LATEX tag is set to YES, doxygen will generate the necessary
+# Makefile rules, Perl scripts and LaTeX code to be able to generate PDF and DVI
+# output from the Perl module output.
+# The default value is: NO.
+# This tag requires that the tag GENERATE_PERLMOD is set to YES.
+
+PERLMOD_LATEX = NO
+
+# If the PERLMOD_PRETTY tag is set to YES, the Perl module output will be nicely
+# formatted so it can be parsed by a human reader. This is useful if you want to
+# understand what is going on. On the other hand, if this tag is set to NO, the
+# size of the Perl module output will be much smaller and Perl will parse it
+# just the same.
+# The default value is: YES.
+# This tag requires that the tag GENERATE_PERLMOD is set to YES.
+
+PERLMOD_PRETTY = YES
+
+# The names of the make variables in the generated doxyrules.make file are
+# prefixed with the string contained in PERLMOD_MAKEVAR_PREFIX. This is useful
+# so different doxyrules.make files included by the same Makefile don't
+# overwrite each other's variables.
+# This tag requires that the tag GENERATE_PERLMOD is set to YES.
+
+PERLMOD_MAKEVAR_PREFIX =
+
+#---------------------------------------------------------------------------
+# Configuration options related to the preprocessor
+#---------------------------------------------------------------------------
+
+# If the ENABLE_PREPROCESSING tag is set to YES, doxygen will evaluate all
+# C-preprocessor directives found in the sources and include files.
+# The default value is: YES.
+
+ENABLE_PREPROCESSING = YES
+
+# If the MACRO_EXPANSION tag is set to YES, doxygen will expand all macro names
+# in the source code. If set to NO, only conditional compilation will be
+# performed. Macro expansion can be done in a controlled way by setting
+# EXPAND_ONLY_PREDEF to YES.
+# The default value is: NO.
+# This tag requires that the tag ENABLE_PREPROCESSING is set to YES.
+
+MACRO_EXPANSION = NO
+
+# If the EXPAND_ONLY_PREDEF and MACRO_EXPANSION tags are both set to YES then
+# the macro expansion is limited to the macros specified with the PREDEFINED and
+# EXPAND_AS_DEFINED tags.
+# The default value is: NO.
+# This tag requires that the tag ENABLE_PREPROCESSING is set to YES.
+
+EXPAND_ONLY_PREDEF = NO
+
+# If the SEARCH_INCLUDES tag is set to YES, the include files in the
+# INCLUDE_PATH will be searched if a #include is found.
+# The default value is: YES.
+# This tag requires that the tag ENABLE_PREPROCESSING is set to YES.
+
+SEARCH_INCLUDES = YES
+
+# The INCLUDE_PATH tag can be used to specify one or more directories that
+# contain include files that are not input files but should be processed by the
+# preprocessor.
+# This tag requires that the tag SEARCH_INCLUDES is set to YES.
+
+INCLUDE_PATH =
+
+# You can use the INCLUDE_FILE_PATTERNS tag to specify one or more wildcard
+# patterns (like *.h and *.hpp) to filter out the header-files in the
+# directories. If left blank, the patterns specified with FILE_PATTERNS will be
+# used.
+# This tag requires that the tag ENABLE_PREPROCESSING is set to YES.
+
+INCLUDE_FILE_PATTERNS =
+
+# The PREDEFINED tag can be used to specify one or more macro names that are
+# defined before the preprocessor is started (similar to the -D option of e.g.
+# gcc). The argument of the tag is a list of macros of the form: name or
+# name=definition (no spaces). If the definition and the "=" are omitted, "=1"
+# is assumed. To prevent a macro definition from being undefined via #undef or
+# recursively expanded use the := operator instead of the = operator.
+# This tag requires that the tag ENABLE_PREPROCESSING is set to YES.
+
+PREDEFINED =
+
+# If the MACRO_EXPANSION and EXPAND_ONLY_PREDEF tags are set to YES then this
+# tag can be used to specify a list of macro names that should be expanded. The
+# macro definition that is found in the sources will be used. Use the PREDEFINED
+# tag if you want to use a different macro definition that overrules the
+# definition found in the source code.
+# This tag requires that the tag ENABLE_PREPROCESSING is set to YES.
+
+EXPAND_AS_DEFINED =
+
+# If the SKIP_FUNCTION_MACROS tag is set to YES then doxygen's preprocessor will
+# remove all references to function-like macros that are alone on a line, have
+# an all uppercase name, and do not end with a semicolon. Such function macros
+# are typically used for boiler-plate code, and will confuse the parser if not
+# removed.
+# The default value is: YES.
+# This tag requires that the tag ENABLE_PREPROCESSING is set to YES.
+
+SKIP_FUNCTION_MACROS = YES
+
+#---------------------------------------------------------------------------
+# Configuration options related to external references
+#---------------------------------------------------------------------------
+
+# The TAGFILES tag can be used to specify one or more tag files. For each tag
+# file the location of the external documentation should be added. The format of
+# a tag file without this location is as follows:
+# TAGFILES = file1 file2 ...
+# Adding location for the tag files is done as follows:
+# TAGFILES = file1=loc1 "file2 = loc2" ...
+# where loc1 and loc2 can be relative or absolute paths or URLs. See the
+# section "Linking to external documentation" for more information about the use
+# of tag files.
+# Note: Each tag file must have a unique name (where the name does NOT include
+# the path). If a tag file is not located in the directory in which doxygen is
+# run, you must also specify the path to the tagfile here.
+
+TAGFILES =
+
+# When a file name is specified after GENERATE_TAGFILE, doxygen will create a
+# tag file that is based on the input files it reads. See section "Linking to
+# external documentation" for more information about the usage of tag files.
+
+GENERATE_TAGFILE =
+
+# If the ALLEXTERNALS tag is set to YES, all external class will be listed in
+# the class index. If set to NO, only the inherited external classes will be
+# listed.
+# The default value is: NO.
+
+ALLEXTERNALS = NO
+
+# If the EXTERNAL_GROUPS tag is set to YES, all external groups will be listed
+# in the modules index. If set to NO, only the current project's groups will be
+# listed.
+# The default value is: YES.
+
+EXTERNAL_GROUPS = YES
+
+# If the EXTERNAL_PAGES tag is set to YES, all external pages will be listed in
+# the related pages index. If set to NO, only the current project's pages will
+# be listed.
+# The default value is: YES.
+
+EXTERNAL_PAGES = YES
+
+# The PERL_PATH should be the absolute path and name of the perl script
+# interpreter (i.e. the result of 'which perl').
+# The default file (with absolute path) is: /usr/bin/perl.
+
+PERL_PATH = /usr/bin/perl
+
+#---------------------------------------------------------------------------
+# Configuration options related to the dot tool
+#---------------------------------------------------------------------------
+
+# If the CLASS_DIAGRAMS tag is set to YES, doxygen will generate a class diagram
+# (in HTML and LaTeX) for classes with base or super classes. Setting the tag to
+# NO turns the diagrams off. Note that this option also works with HAVE_DOT
+# disabled, but it is recommended to install and use dot, since it yields more
+# powerful graphs.
+# The default value is: YES.
+
+CLASS_DIAGRAMS = YES
+
+# You can define message sequence charts within doxygen comments using the \msc
+# command. Doxygen will then run the mscgen tool (see:
+# http://www.mcternan.me.uk/mscgen/)) to produce the chart and insert it in the
+# documentation. The MSCGEN_PATH tag allows you to specify the directory where
+# the mscgen tool resides. If left empty the tool is assumed to be found in the
+# default search path.
+
+MSCGEN_PATH =
+
+# You can include diagrams made with dia in doxygen documentation. Doxygen will
+# then run dia to produce the diagram and insert it in the documentation. The
+# DIA_PATH tag allows you to specify the directory where the dia binary resides.
+# If left empty dia is assumed to be found in the default search path.
+
+DIA_PATH =
+
+# If set to YES the inheritance and collaboration graphs will hide inheritance
+# and usage relations if the target is undocumented or is not a class.
+# The default value is: YES.
+
+HIDE_UNDOC_RELATIONS = YES
+
+# If you set the HAVE_DOT tag to YES then doxygen will assume the dot tool is
+# available from the path. This tool is part of Graphviz (see:
+# http://www.graphviz.org/), a graph visualization toolkit from AT&T and Lucent
+# Bell Labs. The other options in this section have no effect if this option is
+# set to NO
+# The default value is: YES.
+
+HAVE_DOT = YES
+
+# The DOT_NUM_THREADS specifies the number of dot invocations doxygen is allowed
+# to run in parallel. When set to 0 doxygen will base this on the number of
+# processors available in the system. You can set it explicitly to a value
+# larger than 0 to get control over the balance between CPU load and processing
+# speed.
+# Minimum value: 0, maximum value: 32, default value: 0.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+DOT_NUM_THREADS = 0
+
+# When you want a differently looking font in the dot files that doxygen
+# generates you can specify the font name using DOT_FONTNAME. You need to make
+# sure dot is able to find the font, which can be done by putting it in a
+# standard location or by setting the DOTFONTPATH environment variable or by
+# setting DOT_FONTPATH to the directory containing the font.
+# The default value is: Helvetica.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+DOT_FONTNAME = Helvetica
+
+# The DOT_FONTSIZE tag can be used to set the size (in points) of the font of
+# dot graphs.
+# Minimum value: 4, maximum value: 24, default value: 10.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+DOT_FONTSIZE = 10
+
+# By default doxygen will tell dot to use the default font as specified with
+# DOT_FONTNAME. If you specify a different font using DOT_FONTNAME you can set
+# the path where dot can find it using this tag.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+DOT_FONTPATH =
+
+# If the CLASS_GRAPH tag is set to YES then doxygen will generate a graph for
+# each documented class showing the direct and indirect inheritance relations.
+# Setting this tag to YES will force the CLASS_DIAGRAMS tag to NO.
+# The default value is: YES.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+CLASS_GRAPH = YES
+
+# If the COLLABORATION_GRAPH tag is set to YES then doxygen will generate a
+# graph for each documented class showing the direct and indirect implementation
+# dependencies (inheritance, containment, and class references variables) of the
+# class with other documented classes.
+# The default value is: YES.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+COLLABORATION_GRAPH = YES
+
+# If the GROUP_GRAPHS tag is set to YES then doxygen will generate a graph for
+# groups, showing the direct groups dependencies.
+# The default value is: YES.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+GROUP_GRAPHS = YES
+
+# If the UML_LOOK tag is set to YES, doxygen will generate inheritance and
+# collaboration diagrams in a style similar to the OMG's Unified Modeling
+# Language.
+# The default value is: NO.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+UML_LOOK = NO
+
+# If the UML_LOOK tag is enabled, the fields and methods are shown inside the
+# class node. If there are many fields or methods and many nodes the graph may
+# become too big to be useful. The UML_LIMIT_NUM_FIELDS threshold limits the
+# number of items for each type to make the size more manageable. Set this to 0
+# for no limit. Note that the threshold may be exceeded by 50% before the limit
+# is enforced. So when you set the threshold to 10, up to 15 fields may appear,
+# but if the number exceeds 15, the total amount of fields shown is limited to
+# 10.
+# Minimum value: 0, maximum value: 100, default value: 10.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+UML_LIMIT_NUM_FIELDS = 10
+
+# If the TEMPLATE_RELATIONS tag is set to YES then the inheritance and
+# collaboration graphs will show the relations between templates and their
+# instances.
+# The default value is: NO.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+TEMPLATE_RELATIONS = NO
+
+# If the INCLUDE_GRAPH, ENABLE_PREPROCESSING and SEARCH_INCLUDES tags are set to
+# YES then doxygen will generate a graph for each documented file showing the
+# direct and indirect include dependencies of the file with other documented
+# files.
+# The default value is: YES.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+INCLUDE_GRAPH = YES
+
+# If the INCLUDED_BY_GRAPH, ENABLE_PREPROCESSING and SEARCH_INCLUDES tags are
+# set to YES then doxygen will generate a graph for each documented file showing
+# the direct and indirect include dependencies of the file with other documented
+# files.
+# The default value is: YES.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+INCLUDED_BY_GRAPH = YES
+
+# If the CALL_GRAPH tag is set to YES then doxygen will generate a call
+# dependency graph for every global function or class method.
+#
+# Note that enabling this option will significantly increase the time of a run.
+# So in most cases it will be better to enable call graphs for selected
+# functions only using the \callgraph command. Disabling a call graph can be
+# accomplished by means of the command \hidecallgraph.
+# The default value is: NO.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+CALL_GRAPH = NO
+
+# If the CALLER_GRAPH tag is set to YES then doxygen will generate a caller
+# dependency graph for every global function or class method.
+#
+# Note that enabling this option will significantly increase the time of a run.
+# So in most cases it will be better to enable caller graphs for selected
+# functions only using the \callergraph command. Disabling a caller graph can be
+# accomplished by means of the command \hidecallergraph.
+# The default value is: NO.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+CALLER_GRAPH = NO
+
+# If the GRAPHICAL_HIERARCHY tag is set to YES then doxygen will graphical
+# hierarchy of all classes instead of a textual one.
+# The default value is: YES.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+GRAPHICAL_HIERARCHY = YES
+
+# If the DIRECTORY_GRAPH tag is set to YES then doxygen will show the
+# dependencies a directory has on other directories in a graphical way. The
+# dependency relations are determined by the #include relations between the
+# files in the directories.
+# The default value is: YES.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+DIRECTORY_GRAPH = YES
+
+# The DOT_IMAGE_FORMAT tag can be used to set the image format of the images
+# generated by dot. For an explanation of the image formats see the section
+# output formats in the documentation of the dot tool (Graphviz (see:
+# http://www.graphviz.org/)).
+# Note: If you choose svg you need to set HTML_FILE_EXTENSION to xhtml in order
+# to make the SVG files visible in IE 9+ (other browsers do not have this
+# requirement).
+# Possible values are: png, png:cairo, png:cairo:cairo, png:cairo:gd, png:gd,
+# png:gd:gd, jpg, jpg:cairo, jpg:cairo:gd, jpg:gd, jpg:gd:gd, gif, gif:cairo,
+# gif:cairo:gd, gif:gd, gif:gd:gd, svg, png:gd, png:gd:gd, png:cairo,
+# png:cairo:gd, png:cairo:cairo, png:cairo:gdiplus, png:gdiplus and
+# png:gdiplus:gdiplus.
+# The default value is: png.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+DOT_IMAGE_FORMAT = png
+
+# If DOT_IMAGE_FORMAT is set to svg, then this option can be set to YES to
+# enable generation of interactive SVG images that allow zooming and panning.
+#
+# Note that this requires a modern browser other than Internet Explorer. Tested
+# and working are Firefox, Chrome, Safari, and Opera.
+# Note: For IE 9+ you need to set HTML_FILE_EXTENSION to xhtml in order to make
+# the SVG files visible. Older versions of IE do not have SVG support.
+# The default value is: NO.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+INTERACTIVE_SVG = NO
+
+# The DOT_PATH tag can be used to specify the path where the dot tool can be
+# found. If left blank, it is assumed the dot tool can be found in the path.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+DOT_PATH =
+
+# The DOTFILE_DIRS tag can be used to specify one or more directories that
+# contain dot files that are included in the documentation (see the \dotfile
+# command).
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+DOTFILE_DIRS =
+
+# The MSCFILE_DIRS tag can be used to specify one or more directories that
+# contain msc files that are included in the documentation (see the \mscfile
+# command).
+
+MSCFILE_DIRS =
+
+# The DIAFILE_DIRS tag can be used to specify one or more directories that
+# contain dia files that are included in the documentation (see the \diafile
+# command).
+
+DIAFILE_DIRS =
+
+# When using plantuml, the PLANTUML_JAR_PATH tag should be used to specify the
+# path where java can find the plantuml.jar file. If left blank, it is assumed
+# PlantUML is not used or called during a preprocessing step. Doxygen will
+# generate a warning when it encounters a \startuml command in this case and
+# will not generate output for the diagram.
+
+PLANTUML_JAR_PATH =
+
+# When using plantuml, the PLANTUML_CFG_FILE tag can be used to specify a
+# configuration file for plantuml.
+
+PLANTUML_CFG_FILE =
+
+# When using plantuml, the specified paths are searched for files specified by
+# the !include statement in a plantuml block.
+
+PLANTUML_INCLUDE_PATH =
+
+# The DOT_GRAPH_MAX_NODES tag can be used to set the maximum number of nodes
+# that will be shown in the graph. If the number of nodes in a graph becomes
+# larger than this value, doxygen will truncate the graph, which is visualized
+# by representing a node as a red box. Note that doxygen if the number of direct
+# children of the root node in a graph is already larger than
+# DOT_GRAPH_MAX_NODES then the graph will not be shown at all. Also note that
+# the size of a graph can be further restricted by MAX_DOT_GRAPH_DEPTH.
+# Minimum value: 0, maximum value: 10000, default value: 50.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+DOT_GRAPH_MAX_NODES = 50
+
+# The MAX_DOT_GRAPH_DEPTH tag can be used to set the maximum depth of the graphs
+# generated by dot. A depth value of 3 means that only nodes reachable from the
+# root by following a path via at most 3 edges will be shown. Nodes that lay
+# further from the root node will be omitted. Note that setting this option to 1
+# or 2 may greatly reduce the computation time needed for large code bases. Also
+# note that the size of a graph can be further restricted by
+# DOT_GRAPH_MAX_NODES. Using a depth of 0 means no depth restriction.
+# Minimum value: 0, maximum value: 1000, default value: 0.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+MAX_DOT_GRAPH_DEPTH = 0
+
+# Set the DOT_TRANSPARENT tag to YES to generate images with a transparent
+# background. This is disabled by default, because dot on Windows does not seem
+# to support this out of the box.
+#
+# Warning: Depending on the platform used, enabling this option may lead to
+# badly anti-aliased labels on the edges of a graph (i.e. they become hard to
+# read).
+# The default value is: NO.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+DOT_TRANSPARENT = NO
+
+# Set the DOT_MULTI_TARGETS tag to YES to allow dot to generate multiple output
+# files in one run (i.e. multiple -o and -T options on the command line). This
+# makes dot run faster, but since only newer versions of dot (>1.8.10) support
+# this, this feature is disabled by default.
+# The default value is: NO.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+DOT_MULTI_TARGETS = NO
+
+# If the GENERATE_LEGEND tag is set to YES doxygen will generate a legend page
+# explaining the meaning of the various boxes and arrows in the dot generated
+# graphs.
+# The default value is: YES.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+GENERATE_LEGEND = YES
+
+# If the DOT_CLEANUP tag is set to YES, doxygen will remove the intermediate dot
+# files that are used to generate the various graphs.
+# The default value is: YES.
+# This tag requires that the tag HAVE_DOT is set to YES.
+
+DOT_CLEANUP = YES
diff --git a/doc/sphinx/.gitignore b/doc/sphinx/.gitignore
new file mode 100644
index 0000000..04bc145
--- /dev/null
+++ b/doc/sphinx/.gitignore
@@ -0,0 +1 @@
+_*
diff --git a/doc/sphinx/HEJFOG.rst b/doc/sphinx/HEJFOG.rst
new file mode 100644
index 0000000..15bdab4
--- /dev/null
+++ b/doc/sphinx/HEJFOG.rst
@@ -0,0 +1,2 @@
+The HEJ Fixed Order Generator
+=============================
diff --git a/doc/sphinx/Makefile b/doc/sphinx/Makefile
new file mode 100644
index 0000000..71050b8
--- /dev/null
+++ b/doc/sphinx/Makefile
@@ -0,0 +1,20 @@
+# Minimal makefile for Sphinx documentation
+#
+
+# You can set these variables from the command line.
+SPHINXOPTS =
+SPHINXBUILD = python -msphinx
+SPHINXPROJ = reversedHEJ
+SOURCEDIR = .
+BUILDDIR = _build
+
+# Put it first so that "make" without argument is like "make help".
+help:
+ @$(SPHINXBUILD) -M help "$(SOURCEDIR)" "$(BUILDDIR)" $(SPHINXOPTS) $(O)
+
+.PHONY: help Makefile
+
+# Catch-all target: route all unknown targets to Sphinx using the new
+# "make mode" option. $(O) is meant as a shortcut for $(SPHINXOPTS).
+%: Makefile
+ @$(SPHINXBUILD) -M $@ "$(SOURCEDIR)" "$(BUILDDIR)" $(SPHINXOPTS) $(O)
\ No newline at end of file
diff --git a/doc/sphinx/analyses.rst b/doc/sphinx/analyses.rst
new file mode 100644
index 0000000..3a30ca7
--- /dev/null
+++ b/doc/sphinx/analyses.rst
@@ -0,0 +1,145 @@
+Writing custom analyses
+=======================
+
+Both reversed HEJ and HEJ FOG can generate HepMC 3 files, so you can
+always run a rivet analysis on these.
+
+If you want to write your own native analysis, the easiest way is to
+create a new .cc file, for example my_analysis.cc in the
+analysis-plugins/src folder.
+
+An analysis is a class that derives from the abstract :code:`Analysis`
+base class provided by reversed HEJ. It has to implement two public
+functions:
+
+* The :code:`pass_cuts` member function return true if and only if the
+ given event passes the analysis cuts
+
+* The :code:`fill` member function adds an event to the analysis, which
+ for example can be used to fill histograms. Reversed HEJ and HEJ FOG
+ will only give you events for which :code:`pass_cuts` has returned
+ true.
+
+Furthermore, there has to be a global make_analysis function that takes
+the analysis parameters in the form of a YAML :code:`Node` and returns
+a :code:`std::unique_ptr` to the Analysis.
+
+The following code creates the simplest conceivable analysis.::
+
+ #include <memory> // for std::unique_ptr
+
+ #include "RHEJ/Analysis.hh"
+
+ class MyAnalysis: public RHEJ::Analysis {
+ public:
+ MyAnalysis(YAML::Node const & /* config */) {}
+
+ void fill(RHEJ::Event const & /* event */) override{
+ };
+
+ bool pass_cuts(RHEJ::Event const & /* event */) override{
+ return true;
+ };
+
+ };
+
+ extern "C"
+ __attribute__((visibility("default")))
+ std::unique_ptr<RHEJ::Analysis> make_analysis(
+ YAML::Node const & config
+ ){
+ return std::make_unique<MyAnalysis>(config);
+ }
+
+Reversed HEJ and HEJ FOG load analyses dynamically at run time. Since
+the loader cannot deal with the name mangling which is usually done by
+compilers, we have to specify :code:`extern "C"`.
+
+The standard compiler settings for analyses ensure that symbols (for
+example your class name) are not visible in the compiled analysis. This
+ensures that there can be no name clashes with the remaining
+code. However, since reversed HEJ and HEJ FOG have to be able to call
+:code:`make_analysis`, we use
+:code:`__attribute__((visibility("default")))` to ensure its visibility.
+
+To build the analysis go to some build directory and run
+
+.. code-block:: Bash
+
+ cmake some/directory/analysis-plugins
+ make my_analysis
+
+This will create the analysis library libmy_analysis.so in the src
+subdirectory. You can then use it in reversed HEJ or HEJ FOG by adding
+
+.. code-block:: YAML
+
+ analysis:
+ plugin: analysis/build/directory/src/libmy_analysis.so
+
+to the .yml configuration file.
+
+As a more interesting example, here is the code for an analysis that
+sums up the total cross section and prints the result to both standard
+output and a file specified in the .yml config with
+
+.. code-block:: YAML
+
+ analysis:
+ plugin: analysis/build/directory/src/libmy_analysis.so
+ output: outfile
+
+To access the configuration at run time, reversed HEJ uses the yaml-cpp
+library; for more details see the `yaml-cpp tutorial
+<https://github.com/jbeder/yaml-cpp/wiki/Tutorial>`_. The analysis code
+itself is::
+
+ #include <memory>
+ #include <iostream>
+ #include <fstream>
+ #include <string>
+ #include <cmath>
+
+ #include "RHEJ/Analysis.hh"
+ #include "RHEJ/Event.hh"
+
+ #include "yaml-cpp/yaml.h"
+
+ class MyAnalysis: public RHEJ::Analysis {
+ public:
+ MyAnalysis(YAML::Node const & config):
+ xsection_{0.}, xsection_error_{0.},
+ outfile_{config["output"].as<std::string>()}
+ {}
+
+ void fill(RHEJ::Event const & event) override{
+ const double wt = event.central().weight;
+ xsection_ += wt;
+ xsection_error_ += wt*wt;
+ };
+
+ bool pass_cuts(RHEJ::Event const & /* event */) override{
+ return true;
+ };
+
+ ~MyAnalysis(){
+ std::cout << "cross section: " << xsection_ << " +- "
+ << std::sqrt(xsection_error_) << "\n";
+ std::ofstream fout{outfile_};
+ fout << "cross section: " << xsection_ << " +- "
+ << std::sqrt(xsection_error_) << "\n";
+ }
+
+ private:
+ double xsection_, xsection_error_;
+ std::string outfile_;
+
+ };
+
+ extern "C"
+ __attribute__((visibility("default")))
+ std::unique_ptr<RHEJ::Analysis> make_analysis(
+ YAML::Node const & config
+ ){
+ return std::make_unique<MyAnalysis>(config);
+ }
diff --git a/doc/sphinx/conf.py b/doc/sphinx/conf.py
new file mode 100644
index 0000000..4d61cda
--- /dev/null
+++ b/doc/sphinx/conf.py
@@ -0,0 +1,171 @@
+# -*- coding: utf-8 -*-
+#
+# reversed HEJ documentation build configuration file, created by
+# sphinx-quickstart on Fri Sep 15 16:13:57 2017.
+#
+# This file is execfile()d with the current directory set to its
+# containing dir.
+#
+# Note that not all possible configuration values are present in this
+# autogenerated file.
+#
+# All configuration values have a default; values that are commented out
+# serve to show the default.
+
+# If extensions (or modules to document with autodoc) are in another directory,
+# add these directories to sys.path here. If the directory is relative to the
+# documentation root, use os.path.abspath to make it absolute, like shown here.
+#
+# import os
+# import sys
+# sys.path.insert(0, os.path.abspath('.'))
+
+
+# -- General configuration ------------------------------------------------
+
+# If your documentation needs a minimal Sphinx version, state it here.
+#
+# needs_sphinx = '1.0'
+
+# Add any Sphinx extension module names here, as strings. They can be
+# extensions coming with Sphinx (named 'sphinx.ext.*') or your custom
+# ones.
+extensions = ['sphinx.ext.mathjax',
+ 'sphinx.ext.githubpages']
+
+# Add any paths that contain templates here, relative to this directory.
+templates_path = ['_templates']
+
+# The suffix(es) of source filenames.
+# You can specify multiple suffix as a list of string:
+#
+# source_suffix = ['.rst', '.md']
+source_suffix = '.rst'
+
+# The master toctree document.
+master_doc = 'index'
+
+# General information about the project.
+project = u'reversed HEJ'
+copyright = u'2017, Jeppe Andersen, Tuomas Hapola, Andreas Maier, Jennifer Smillie'
+author = u'Jeppe Andersen, Tuomas Hapola, Andreas Maier, Jennifer Smillie'
+
+# The version info for the project you're documenting, acts as replacement for
+# |version| and |release|, also used in various other places throughout the
+# built documents.
+#
+# The short X.Y version.
+version = u'0.0.1'
+# The full version, including alpha/beta/rc tags.
+release = u'0.0.1'
+
+# The language for content autogenerated by Sphinx. Refer to documentation
+# for a list of supported languages.
+#
+# This is also used if you do content translation via gettext catalogs.
+# Usually you set "language" from the command line for these cases.
+language = None
+
+# List of patterns, relative to source directory, that match files and
+# directories to ignore when looking for source files.
+# This patterns also effect to html_static_path and html_extra_path
+exclude_patterns = ['_build', 'Thumbs.db', '.DS_Store']
+
+highlight_language = 'C++'
+
+# The name of the Pygments (syntax highlighting) style to use.
+pygments_style = 'sphinx'
+
+# If true, `todo` and `todoList` produce output, else they produce nothing.
+todo_include_todos = False
+
+
+# -- Options for HTML output ----------------------------------------------
+
+# The theme to use for HTML and HTML Help pages. See the documentation for
+# a list of builtin themes.
+#
+html_theme = 'alabaster'
+
+# Theme options are theme-specific and customize the look and feel of a theme
+# further. For a list of options available for each theme, see the
+# documentation.
+#
+# html_theme_options = {}
+
+# Add any paths that contain custom static files (such as style sheets) here,
+# relative to this directory. They are copied after the builtin static files,
+# so a file named "default.css" will overwrite the builtin "default.css".
+html_static_path = ['_static']
+
+# Custom sidebar templates, must be a dictionary that maps document names
+# to template names.
+#
+# This is required for the alabaster theme
+# refs: http://alabaster.readthedocs.io/en/latest/installation.html#sidebars
+html_sidebars = {
+ '**': [
+ 'about.html',
+ 'navigation.html',
+ 'relations.html', # needs 'show_related': True theme option to display
+ 'searchbox.html',
+ 'donate.html',
+ ]
+}
+
+
+# -- Options for HTMLHelp output ------------------------------------------
+
+# Output file base name for HTML help builder.
+htmlhelp_basename = 'reversedHEJdoc'
+
+
+# -- Options for LaTeX output ---------------------------------------------
+
+latex_elements = {
+ # The paper size ('letterpaper' or 'a4paper').
+ #
+ # 'papersize': 'letterpaper',
+
+ # The font size ('10pt', '11pt' or '12pt').
+ #
+ # 'pointsize': '10pt',
+
+ # Additional stuff for the LaTeX preamble.
+ #
+ # 'preamble': '',
+
+ # Latex figure (float) alignment
+ #
+ # 'figure_align': 'htbp',
+}
+
+# Grouping the document tree into LaTeX files. List of tuples
+# (source start file, target name, title,
+# author, documentclass [howto, manual, or own class]).
+latex_documents = [
+ (master_doc, 'reversedHEJ.tex', u'reversed HEJ Documentation',
+ u'Jeppe Andersen, Tuomas Hapola, Andreas Maier, Jennifer Smillie', 'manual'),
+]
+
+
+# -- Options for manual page output ---------------------------------------
+
+# One entry per manual page. List of tuples
+# (source start file, name, description, authors, manual section).
+man_pages = [
+ (master_doc, 'reversedhej', u'reversed HEJ Documentation',
+ [author], 1)
+]
+
+
+# -- Options for Texinfo output -------------------------------------------
+
+# Grouping the document tree into Texinfo files. List of tuples
+# (source start file, target name, title, author,
+# dir menu entry, description, category)
+texinfo_documents = [
+ (master_doc, 'reversedHEJ', u'reversed HEJ Documentation',
+ author, 'reversedHEJ', 'One line description of project.',
+ 'Miscellaneous'),
+]
diff --git a/doc/sphinx/index.rst b/doc/sphinx/index.rst
new file mode 100644
index 0000000..d00f701
--- /dev/null
+++ b/doc/sphinx/index.rst
@@ -0,0 +1,23 @@
+.. reversed HEJ documentation master file, created by
+ sphinx-quickstart on Fri Sep 15 16:13:57 2017.
+
+Welcome to reversed HEJ's documentation!
+========================================
+
+.. toctree::
+ :maxdepth: 2
+ :caption: Contents:
+
+ intro
+ installation
+ rHEJ
+ HEJFOG
+ analyses
+
+
+Indices and tables
+==================
+
+* :ref:`genindex`
+* :ref:`modindex`
+* :ref:`search`
diff --git a/doc/sphinx/installation.rst b/doc/sphinx/installation.rst
new file mode 100644
index 0000000..11e4437
--- /dev/null
+++ b/doc/sphinx/installation.rst
@@ -0,0 +1,2 @@
+Installation
+============
diff --git a/doc/sphinx/intro.rst b/doc/sphinx/intro.rst
new file mode 100644
index 0000000..c516b33
--- /dev/null
+++ b/doc/sphinx/intro.rst
@@ -0,0 +1,2 @@
+Introduction
+============
diff --git a/doc/sphinx/rHEJ.rst b/doc/sphinx/rHEJ.rst
new file mode 100644
index 0000000..36934b4
--- /dev/null
+++ b/doc/sphinx/rHEJ.rst
@@ -0,0 +1,2 @@
+Running reversed HEJ
+====================
diff --git a/include/RHEJ/Analysis.hh b/include/RHEJ/Analysis.hh
index 04d7eef..06dc12b 100644
--- a/include/RHEJ/Analysis.hh
+++ b/include/RHEJ/Analysis.hh
@@ -1,86 +1,91 @@
/** \file Analysis.hh
* \brief Header file for the Analysis Interface (ABC)
*
* Contains the Analysis struct and the Dynamic Analysis
- * struct. Also contains some histogram centric functions
+ * struct. Also contains some histogram centric functions
* such as histogram_names(), parameters(), draw_histogram()
* update_parameter()
*/
#pragma once
#include <vector>
#include <string>
+namespace YAML{
+ class Node;
+}
+
//! Reversed HEJ Namespace
namespace RHEJ{
class Event;
/** \struct Analysis Analysis.hh "include/RHEJ/Analysis.hh"
* \brief Analysis Struct
- *
- * Contains a fill function which needs to be overridden by any
- * inheriting class before it can be utilized.
+ *
+ * This is the interface that all analyses should implement,
+ * i.e. all custom analyses have to derived from this struct.
*/
struct Analysis{
+ //! Fill event into analysis (e.g. to histograms)
virtual void fill(Event const &) = 0;
+ //! Decide whether an event passes the cuts
+ virtual bool pass_cuts(Event const &) = 0;
virtual ~Analysis() = default;
};
/** \struct DynamicAnalysis Analysis.hh "include/RHEJ/Analysis.hh"
* \brief DynamicAnalysis struct to be used with HEJview
*
* Inherits from the standard Analysis struct.
*/
struct DynamicAnalysis : Analysis {
/** \struct Parameter Analysis.hh "include/RHEJ/Analysis.hh"
* \brief Parameter Struct which contains Parameters of Analysis
*
* Contains all info necessary about a parameter of the analysis.
*/
struct Parameter{
std::string name; /**< Name of parameter */
double value; /**< Value of Paramter */
double min; /**< Min Value allowed for parameter */
double max; /**< Max Value allowed for parameter */
};
//! Virtual Histogram Names Function
- /**
+ /**
* Returns a vector of the names of the histograms.
*/
virtual std::vector<std::string> histogram_names() const = 0;
//! Virtual Parameters Function
- /**
- * Returns a vector which contains the parameters to be placed into the
- * histograms
+ /**
+ * Returns the a posteriori adjustable parameters of the analysis,
+ * for instance weight-weight cuts
*/
virtual std::vector<Parameter> parameters() const = 0;
//! Virtual Draw Histogram Function
- /**
+ /**
* @param name Name of Histogram to be drawn
* @param filename Name of Output file containing Histogram
*
- * This function will update the appropriate histogram?
*/
virtual void draw_histogram(
std::string const & name, std::string const & filename
) const = 0;
//! Virtual Update Parameter Function
/**
* @param name Name of the Parameter which is to be updated
* @param value Value to be updated to
- * @returns 0 or 1 depending on if value is updated?
+ * @returns true if any histograms have changed
*
- * How is its decided whether the value will be updated?
*/
virtual bool update_parameter(std::string const & name, double value) = 0;
//! Dynamic Analysis Destructor
virtual ~DynamicAnalysis() = default;
};
}
diff --git a/include/RHEJ/EmptyAnalysis.hh b/include/RHEJ/EmptyAnalysis.hh
index 1733892..2c879ba 100644
--- a/include/RHEJ/EmptyAnalysis.hh
+++ b/include/RHEJ/EmptyAnalysis.hh
@@ -1,31 +1,32 @@
/** \file EmptyAnalysis.hh
* \brief Details the EmptyAnalysis Class. (The trivial Analysis case.)
*/
#pragma once
#include <memory>
#include "RHEJ/Analysis.hh"
//! YAML Namespace
namespace YAML{
class Node;
}
namespace RHEJ{
/** \struct EmptyAnalysis EmptyAnalysis.hh "include/RHEJ/EmptyAnalysis.hh"
* \brief EmptyAnalysis struct. Trivial Analysis case.
*
* Inherits from the Analysis Struct and defines its own fill function (as necessary)
* This function does nothing (which is the trivial case.
*/
struct EmptyAnalysis: Analysis{
static std::unique_ptr<Analysis> create(YAML::Node const & parameters);
virtual void fill(Event const &) override;
+ virtual bool pass_cuts(Event const &) override;
virtual ~EmptyAnalysis() override = default;
};
}
diff --git a/include/RHEJ/MatrixElement.hh b/include/RHEJ/MatrixElement.hh
index fac8e30..8da4a3e 100644
--- a/include/RHEJ/MatrixElement.hh
+++ b/include/RHEJ/MatrixElement.hh
@@ -1,219 +1,226 @@
/** \file MatrixElement.hh
* \brief The header file which contains the MatrixElement Class
*
* This contains the MatrixElement Class which contains many functions
* used to calculate many different MatrixElements and their components.
*/
#pragma once
+#include <functional>
+
#include "RHEJ/utility.hh"
#include "RHEJ/HiggsCouplingSettings.hh"
#include "CLHEP/Vector/LorentzVector.h"
namespace RHEJ{
/** \class MatrixElement MatrixElement.hh "include/RHEJ/MatrixElement.hh
* \brief MatrixElement class. Functions for obtaining various ME and components.
*/
class MatrixElement{
public:
/** \brief MatrixElement Constructor
+ * @param alpha_s Function taking the renormalisation scale
+ * and returning the strong coupling constant
+ * @param jet_def Jet definition for deciding which particles
+ * are members of extremal jets
+ * @param log_corr Whether to including logarithmic corrections
+ * from the running of the coupling
+ * @param Hgg_settings Settings for the Higgs coupling to gluons
*/
MatrixElement(
+ std::function<double (double)> alpha_s,
fastjet::JetDefinition jet_def, double jetptmin,
bool log_corr,
HiggsCouplingSettings Hgg_settings
);
/**
* \brief regulated HEJ matrix element
- * @param alpha_s Value of the strong coupling constant
* @param mur Value of the renormalisation scale
* @param incoming Incoming particles
* @param partons Outgoing particles
* @param check_momenta Special treatment for partons inside extremal jets
* @returns The HEJ matrix element including virtual corrections
*
* cf. eq. (22) in \ref Andersen:2011hs
* Incoming particles should be ordered by ascending z momentum.
* Outgoing particles should be ordered by ascending rapidity.
*
* \internal Relation to standard HEJ Met2: MatrixElement = Met2*shat^2/(pdfta*pdftb)
*/
double operator()(
- double alpha_s, double mur,
+ double mur,
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> const & outgoing,
bool check_momenta
) const;
//! HEJ tree-level matrix element
/**
- * @param alpha_s Value of the strong coupling constant
* @param mur Value of the renormalisation scale
* @param incoming Incoming particles
* @param outgoing Outgoing particles
* @param check_momenta Special treatment for partons inside extremal jets
* @returns The HEJ matrix element without virtual corrections
*
* cf. eq. (22) in \ref Andersen:2011hs
* Incoming particles should be ordered by ascending z momentum.
* Outgoing particles should be ordered by ascending rapidity.
*/
double tree(
- double alpha_s, double mur,
+ double mur,
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> const & outgoing,
bool check_momenta
) const;
//! HEJ tree-level matrix element - parametric part
/**
- * @param alpha_s Value of the strong coupling constant
* @param mur Value of the renormalisation scale
* @param outgoing Outgoing particles
* @returns The parametric part of the tree matrix element
*
* cf. eq. (22) in \ref Andersen:2011hs
*
* The tree level matrix element factorises into a parametric part
* which depends on the theory parameters (alpha_s and scale)
* and a kinematic part comprising the dependence on the particle momenta
* and colour factors. This function returns the former.
*/
double tree_param(
- double alpha_s, double mur,
+ double mur,
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> const & outgoing
) const;
//! HEJ tree-level matrix element - kinematic part
/**
* @param incoming Incoming particles
* @param partons Outgoing particles
* @param check_momenta Special treatment for partons inside extremal jets
* @returns The kinematic part of the tree matrix element
*
* cf. eq. (22) in \ref Andersen:2011hs
* Incoming particles should be ordered by ascending z momentum.
* Outgoing particles should be ordered by ascending rapidity.
*
* The tree level matrix element factorises into a parametric part
* which depends on the theory parameters (alpha_s and scale)
* and a kinematic part comprising the dependence on the particle momenta
* and colour factors. This function returns the latter.
*/
double tree_kin(
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> const & outgoing,
bool check_momenta
) const;
/**
* \brief Calculates the Virtual Corrections
- * @param alpha_s Value of the strong coupling constant
* @param mur Value of the renormalisation scale
* @param in Incoming particles
* @param out Outgoing particles
* @returns The Virtual Corrections of the Matrix Element
*
* Incoming particles should be ordered by ascending z momentum.
* Outgoing particles should be ordered by ascending rapidity.
*
* The all order virtual corrections to LL in the MRK limit is
* given by replacing 1/t in the scattering amplitude according to the
* lipatov ansatz.
*
* cf. second-to-last line of eq. (22) in \ref Andersen:2011hs
* note that indices are off by one, i.e. out[0].p corresponds to p_1
*/
double virtual_corrections(
- double alpha_s, double mur,
+ double mur,
std::array<Sparticle, 2> const & in,
std::vector<Sparticle> const & out
) const;
private:
/**
* \brief cf. last line of eq. (22) in \ref Andersen:2011hs
- * @param alpha_s Value of Strong Coupling Constant
* @param mur Value of Renormalisation Scale
* @param q_j ???
* @param lambda ???
*/
double omega0(
double alpha_s, double mur,
fastjet::PseudoJet const & q_j, double lambda
) const;
double tree_kin_jets(
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> partons,
bool check_momenta
) const;
double tree_kin_Higgs(
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> const & outgoing,
bool check_momenta
) const;
double tree_kin_Higgs_first(
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> const & outgoing,
bool check_momenta
) const;
double tree_kin_Higgs_last(
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> const & outgoing,
bool check_momenta
) const;
/**
* \brief Higgs inbetween extremal partons.
*
* Note that in the case of unordered emission, the Higgs is *always*
* treated as if in between the extremal (FKL) partons, even if its
* rapidity is outside the extremal parton rapidities
*/
double tree_kin_Higgs_between(
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> const & outgoing,
bool check_momenta
) const;
double tree_param_partons(
double alpha_s, double mur,
std::vector<Sparticle> const & partons
) const;
std::vector<int> in_extremal_jet_indices(
std::vector<fastjet::PseudoJet> const & partons
) const;
std::vector<Sparticle> tag_extremal_jet_partons(
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> out_partons, bool check_momenta
) const;
double MH2_forwardH(
RHEJ::ParticleID id,
CLHEP::HepLorentzVector p1out, CLHEP::HepLorentzVector p1in,
CLHEP::HepLorentzVector p2out, CLHEP::HepLorentzVector p2in,
CLHEP::HepLorentzVector pH,
double t1, double t2
) const;
+ std::function<double (double)> alpha_s_;
+
HiggsCouplingSettings Hgg_settings_;
fastjet::JetDefinition jet_def_;
double jetptmin_;
bool log_corr_;
};
}
diff --git a/include/RHEJ/config.hh b/include/RHEJ/config.hh
index bab25fb..5295d5c 100644
--- a/include/RHEJ/config.hh
+++ b/include/RHEJ/config.hh
@@ -1,114 +1,343 @@
/** \file config.hh
* \brief The file which handles the configuration file parameters
*
* Contains the JetParameters Struct and ScaleConfig Struct. Also
* contains the ScaleGenerator Class and the EventTreatment class.
* Config struct is also defined here. The configuration files
* parameters are read and then stored within this objects.
*/
#pragma once
#include <string>
#include <memory>
#include "fastjet/JetDefinition.hh"
#include "RHEJ/ScaleFunction.hh"
#include "RHEJ/optional.hh"
#include "RHEJ/output_formats.hh"
#include "RHEJ/event_classes.hh"
#include "RHEJ/HiggsCouplingSettings.hh"
+#include "RHEJ/exceptions.hh"
+#include "RHEJ/utility.hh"
#include "yaml-cpp/yaml.h"
namespace RHEJ{
/**! \struct JetParameters config.hh "include/RHEJ/config.hh"
* \brief Contains Jet Definition and Minimum Pt to be a Jet.
*
* Contains the Fastjet definition of a Jet and a threshold (minimum) value for pt
*/
struct JetParameters{
fastjet::JetDefinition def; /**< Jet Definition */
double min_pt; /**< Minimum Jet Pt */
};
/**! \struct ScaleConfig config.hh "include/RHEJ/config.hh"
* \brief Sets up the possible choices for scale variation?
*
* There is a vector which contains the possible scale choices and
* also vector of doubles with the scale factors along with a max
* scale ratio.
*/
struct ScaleConfig{
std::vector<std::unique_ptr<ScaleFunction>> scales; /**< Vector of possible Scale choices */
std::vector<double> scale_factors; /**< Vector of possible Scale Factors */
double max_scale_ratio; /**< Maximum ratio for the scales */
};
/**! \class ScaleGenerator config.hh "include/RHEJ/config.hh"
* \brief A Class used to generate scales
*
* This class has a clustered event class and scale config struct
* within it.
*/
class ScaleGenerator{
public:
ScaleGenerator() = default;
explicit ScaleGenerator(ScaleConfig && settings):
settings_{std::move(settings)}
{}
Event operator()(Event ev) const;
ScaleConfig const & settings() const;
private:
ScaleConfig settings_;
};
/**
* \brief Enumeration which deals with how to treat events.
*
* The program will decide on how to treat an event based on
* the value of this enumeration.
*/
enum class EventTreatment{
reweight, /**< Reweigh the event */
keep, /**< Keep the event */
discard, /**< Discard the event */
};
/**
* \brief An ordered Map with EventClass keys and mapped values EventTreatment
*
* If called upon with EventTreatMap.at(EventClass) it will return the treatment which has
* been specified for that EventClass.
*/
using EventTreatMap = std::map<event_class::EventClass, EventTreatment>;
/**! \struct Config config.hh "include/RHEJ/config.hh"
* \brief Config Struct for user input parameters.
*
* The struct which handles the input card given by the user.
+ *
+ * \internal To add a new option:
+ * 1. Add a member to the Config struct.
+ * 2. Add the option name to the "supported" Node in
+ * get_supported_options.
+ * 3. Initialise the new Config member in to_Config.
+ * The functions set_from_yaml (for mandatory options) and
+ * set_from_yaml_if_defined (non-mandatory) may be helpful.
+ * 4. Add a new entry (with short description) to config.yaml
*/
struct Config {
ScaleGenerator scale_gen; /**< Scale */
JetParameters resummation_jets; /**< Resummation Jets */
JetParameters fixed_order_jets; /**< Fixed-Order Jets */
double min_extparton_pt; /**< Min External Parton Pt*/
double max_ext_soft_pt_fraction; /**< Min External Parton Pt Fraction */
int trials; /**< Number of resummation points to generate per FO */
bool log_correction; /**< Log Correct On or Off */
bool unweight; /**< Unweight On or Off */
std::vector<OutputFile> output; /**< Output File Type */
optional<std::string> RanLux_init; /**< Ranlux File Choice*/
EventTreatMap treat; /**< Event Treat Map? */
YAML::Node analysis_parameters; /**< Analysis Parameters */
//! Settings for Effective Higgs-Gluon Coupling
HiggsCouplingSettings Higgs_coupling;
};
+ //! Load configuration from file
+ /**
+ * @param config_file Name of the YAML configuration file
+ * @returns The reversed HEJ 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
+ );
+
+ JetParameters get_jet_parameters(
+ YAML::Node const & node, std::string const & entry
+ );
+
+ HiggsCouplingSettings get_Higgs_coupling(
+ YAML::Node const & node, std::string const & entry
+ );
+
+ ScaleConfig to_ScaleConfig(YAML::Node const & yaml);
+
+ 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);
+
+ 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...) + ": value " + ex.what()
+ + " cannot be converted to a " + type_string(setting)
+ };
+ }
+ 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<RHEJ::OutputFile> {
+ static Node encode(RHEJ::OutputFile const & outfile);
+ static bool decode(Node const & node, RHEJ::OutputFile & out);
+ };
}
diff --git a/include/RHEJ/utility.hh b/include/RHEJ/utility.hh
index 542e710..9caab71 100644
--- a/include/RHEJ/utility.hh
+++ b/include/RHEJ/utility.hh
@@ -1,118 +1,141 @@
/**
* \file utility.hh
* Author: Tuomas Hapola
* \brief Contains the struct Sparticle which contains particle information and parameters.
*
* Also contains some functions which are useful for comparison on particles parameters.
*/
#pragma once
#include <algorithm>
+#include <boost/core/demangle.hpp>
+
// FastJet Includes
#include "fastjet/PseudoJet.hh"
#include "fastjet/ClusterSequence.hh"
#include "RHEJ/PDG_codes.hh"
namespace RHEJ{
/** \struct Sparticle utility.hh "include/RHEJ/utility.hh"
* \brief The struct which contains all the parameters of a particle.
*
* Constituted by PiD and PseudoJet Objects.
*/
struct Sparticle {
ParticleID type;
fastjet::PseudoJet p;
//! get rapidity function
double rapidity() const{
return p.rapidity();
}
//! get perp function
double perp() const{
return p.perp();
}
//! get Px function
double px() const{
return p.px();
}
//! get Py function
double py() const{
return p.py();
}
//! Get Pz function
double pz() const{
return p.pz();
}
//! Get P0 function
double E() const{
return p.E();
}
};
/** \struct rapidity_less utility.hh "include/RHEJ/utility.hh
- * \brief A struct which allows quick comparison of two different jets rapidities.
+ * \brief A struct which allows quick comparison of two different jets rapidities.
*/
struct rapidity_less{
template<class FourVector>
bool operator()(FourVector const & p1, FourVector const & p2){
return p1.rapidity() < p2.rapidity();
}
};
-
+
/** \struct pz_less utility.hh "include/RHEJ/utility.hh
- * \brief A struct which allows quick comparison of Pz between two jets.
+ * \brief A struct which allows quick comparison of Pz between two jets.
*/
struct pz_less{
template<class FourVector>
bool operator()(FourVector const & p1, FourVector const & p2){
return p1.pz() < p2.pz();
}
};
-
+
inline
std::vector<fastjet::PseudoJet> to_PseudoJet(
std::vector<Sparticle> const & v
){
std::vector<fastjet::PseudoJet> result;
for(auto && sp: v) result.emplace_back(sp.p);
return result;
}
inline
bool is_parton(Sparticle const & p){
return is_parton(p.type);
}
inline bool is_AWZH_boson(Sparticle const & particle){
return is_AWZH_boson(particle.type);
}
inline
std::vector<Sparticle> filter_partons(
std::vector<Sparticle> const & v
){
std::vector<Sparticle> result;
result.reserve(v.size());
std::copy_if(
begin(v), end(v), std::back_inserter(result),
[](Sparticle const & p){ return is_parton(p); }
);
return result;
}
template<class T, class... Args>
std::unique_ptr<T> make_unique(Args&&... a){
return std::unique_ptr<T>{new T{std::forward<Args>(a)...}};
}
template<typename T, typename... U>
constexpr
std::array<T, 1 + sizeof...(U)> make_array(T t, U&&... rest){
return {{t, std::forward<U>(rest)...}};
}
+ inline
+ std::string join(
+ std::string const & /* delim */, std::string const & str
+ ){
+ return str;
+ }
+
+ template<typename... Strings>
+ std::string join(
+ std::string const & delim,
+ std::string const & first, std::string const & second,
+ Strings&&... rest
+ ){
+ return join(delim, first + delim + second, std::forward<Strings>(rest)...);
+ }
+
+ template<typename T>
+ std::string type_string(T&&){
+ return boost::core::demangle(typeid(T).name());
+ }
+
}
diff --git a/src/EmptyAnalysis.cc b/src/EmptyAnalysis.cc
index 0ebf987..3eb579e 100644
--- a/src/EmptyAnalysis.cc
+++ b/src/EmptyAnalysis.cc
@@ -1,54 +1,58 @@
#include "RHEJ/EmptyAnalysis.hh"
#include "RHEJ/exceptions.hh"
#include <iostream>
#include <string>
#include <vector>
#include "yaml-cpp/yaml.h"
namespace RHEJ{
namespace{
std::vector<std::string> param_as_strings(YAML::Node const & parameters){
using YAML::NodeType;
switch(parameters.Type()){
case NodeType::Null:
case NodeType::Undefined:
return {};
case NodeType::Scalar:
return {parameters.as<std::string>()};
case NodeType::Sequence: {
std::vector<std::string> param_strings;
for(auto && param: parameters){
param_strings.emplace_back(param.as<std::string>());
}
return param_strings;
}
case NodeType::Map: {
std::vector<std::string> param_strings;
for(auto && param: parameters){
param_strings.emplace_back(param.first.as<std::string>());
}
return param_strings;
}
default:;
}
throw std::logic_error{"unreachable"};
}
}
std::unique_ptr<Analysis> EmptyAnalysis::create(
YAML::Node const & parameters
){
const auto param_strings = param_as_strings(parameters);
if(! param_strings.empty()){
std::string error{"Unknown analysis parameter(s):"};
for(auto && p: param_strings) error += " " + p;
throw unknown_option{error};
}
return std::unique_ptr<Analysis>{new EmptyAnalysis{}};
}
void EmptyAnalysis::fill(Event const &){
// do nothing
}
+
+ bool EmptyAnalysis::pass_cuts(Event const &){
+ return true;
+ }
}
diff --git a/src/EventReweighter.cc b/src/EventReweighter.cc
index 717524b..8bd2c2f 100644
--- a/src/EventReweighter.cc
+++ b/src/EventReweighter.cc
@@ -1,320 +1,323 @@
#include "RHEJ/EventReweighter.hh"
#include <string>
#include <unordered_map>
#include "RHEJ/PhaseSpacePoint.hh"
#include "RHEJ/PDG_codes.hh"
#include "RHEJ/config.hh" // ScaleGenerator definition
#include "RHEJ/debug.hh"
namespace RHEJ{
using EventClass = event_class::EventClass;
namespace {
static_assert(
std::numeric_limits<double>::has_quiet_NaN,
"no quiet NaN for double"
);
constexpr double NaN = std::numeric_limits<double>::quiet_NaN();
UnclusteredEvent to_UnclusteredEvent(PhaseSpacePoint const & psp){
UnclusteredEvent result;
result.incoming = psp.incoming();
std::sort(
begin(result.incoming), end(result.incoming),
[](Sparticle o1, Sparticle o2){return o1.p.pz()<o2.p.pz();}
);
assert(result.incoming.size() == 2);
result.outgoing = psp.outgoing();
assert(
std::is_sorted(
begin(result.outgoing), end(result.outgoing),
rapidity_less{}
)
);
assert(result.outgoing.size() >= 2);
result.central.mur = NaN;
result.central.muf = NaN;
result.central.weight = psp.weight();
return result;
}
} // anonymous namespace
EventReweighter::EventReweighter(
LHEF::HEPRUP const & heprup,
fastjet::JetDefinition jet_def, double jetptmin,
EventTreatMap treat,
double extpartonptmin, double max_ext_soft_pt_fraction,
bool log_corr,
HiggsCouplingSettings Higgs_coupling
):
EventReweighter{
RHEJ::Beam{
heprup.EBMUP.first,
{{
static_cast<RHEJ::ParticleID>(heprup.IDBMUP.first),
static_cast<RHEJ::ParticleID>(heprup.IDBMUP.second)
}}
},
heprup.PDFSUP.first,
jet_def, jetptmin,
std::move(treat),
extpartonptmin, max_ext_soft_pt_fraction,
log_corr,
std::move(Higgs_coupling)
}
{
if(heprup.EBMUP.second != E_beam_){
throw std::invalid_argument(
"asymmetric beam: " + std::to_string(E_beam_)
+ " ---> <--- " + std::to_string(heprup.EBMUP.second)
);
};
if(heprup.PDFSUP.second != pdf_.id()){
throw std::invalid_argument(
"conflicting PDF ids: " + std::to_string(pdf_.id())
+ " vs. " + std::to_string(heprup.PDFSUP.second)
);
}
}
EventReweighter::EventReweighter(
Beam beam,
int pdf_id,
fastjet::JetDefinition jet_def, double jetptmin,
EventTreatMap treat,
double extpartonptmin, double max_ext_soft_pt_fraction,
bool log_corr,
HiggsCouplingSettings Higgs_coupling
):
extpartonptmin_{extpartonptmin},
max_ext_soft_pt_fraction_{max_ext_soft_pt_fraction},
jetptmin_{jetptmin},
E_beam_{beam.E},
jet_def_{jet_def},
pdf_{pdf_id, beam.type.front(), beam.type.back()},
- MEt2_{jet_def, jetptmin, log_corr, std::move(Higgs_coupling)},
+ MEt2_{
+ [this](double mu){ return pdf_.Halphas(mu); },
+ jet_def, jetptmin,
+ log_corr, std::move(Higgs_coupling)
+ },
treat_{std::move(treat)}
{}
PDF const & EventReweighter::pdf() const{
return pdf_;
}
EventClass EventReweighter::last_event_class() const{
return last_event_class_;
}
namespace{
static_assert(
std::numeric_limits<double>::has_infinity, "infinity not supported"
);
// there is no simple way to create a vector of move-only objects
// this is a clumsy work-around
ScaleGenerator make_identity(){
std::vector<std::unique_ptr<ScaleFunction>> scale_fun;
scale_fun.emplace_back(new InputScales);
return ScaleGenerator{
ScaleConfig{
std::move(scale_fun), {}, std::numeric_limits<double>::infinity()
}
};
}
const ScaleGenerator identity{make_identity()};
}
std::vector<Event> EventReweighter::reweight(
Event const & input_ev, int num_events
){
return reweight(input_ev, num_events, identity);
}
std::vector<Event> EventReweighter::reweight(
Event const & input_ev, int num_events,
ScaleGenerator const & gen_scales
){
auto res_events = gen_res_events(input_ev, num_events);
if(res_events.empty()) return {};
for(auto & event: res_events) event = gen_scales(event);
return rescale(input_ev, std::move(res_events));
}
/**
* \brief main generation/reweighting function: generate phase space points and divide out Born factors
*/
std::vector<Event> EventReweighter::gen_res_events(
Event const & ev,
int phase_space_points
){
assert(ev.variations().empty());
last_event_class_ = classify(ev);
switch(treat_.at(last_event_class_)){
case EventTreatment::discard: return {};
case EventTreatment::keep:
if(! jets_pass_resummation_cuts(ev)) return {};
else return {ev};
default:;
}
const double Born_shat = shat(ev);
std::vector<Event> resummation_events;
for(int psp_number = 0; psp_number < phase_space_points; ++psp_number){
PhaseSpacePoint psp{
ev,
jet_def_, jetptmin_,
extpartonptmin_, max_ext_soft_pt_fraction_
};
if(psp.weight() == 0.) continue;
resummation_events.emplace_back(
to_UnclusteredEvent(std::move(psp)), jet_def_, jetptmin_
);
auto & new_event = resummation_events.back();
assert(new_event.variations().empty());
new_event.central().mur = ev.central().mur;
new_event.central().muf = ev.central().muf;
const double resum_shat = shat(new_event);
new_event.central().weight *= ev.central().weight*Born_shat*Born_shat/
(phase_space_points*resum_shat*resum_shat);
}
return resummation_events;
}
std::vector<Event> EventReweighter::rescale(
Event const & Born_ev,
std::vector<Event> events
) const{
const double Born_pdf = pdf_factors(Born_ev).central;
const double Born_ME = tree_matrix_element(Born_ev);
for(auto & cur_event: events){
const auto pdf = pdf_factors(cur_event);
assert(pdf.variations.size() == cur_event.variations().size());
const auto ME = matrix_elements(cur_event);
assert(ME.variations.size() == cur_event.variations().size());
cur_event.central().weight *= pdf.central*ME.central/(Born_pdf*Born_ME);
for(size_t i = 0; i < cur_event.variations().size(); ++i){
cur_event.variations(i).weight *=
pdf.variations[i]*ME.variations[i]/(Born_pdf*Born_ME);
}
}
return events;
};
bool EventReweighter::jets_pass_resummation_cuts(
Event const & ev
) const{
const auto out_as_PseudoJet = to_PseudoJet(filter_partons(ev.outgoing()));
fastjet::ClusterSequence cs{out_as_PseudoJet, jet_def_};
return cs.inclusive_jets(jetptmin_).size() == ev.jets().size();
}
EventReweighter::EventFactors
EventReweighter::pdf_factors(Event const & ev) const{
auto const & a = ev.incoming().front();
auto const & b = ev.incoming().back();
const double xa = a.p.e()/E_beam_;
const double xb = b.p.e()/E_beam_;
EventFactors result;
std::unordered_map<double, double> known_pdf;
result.central =
pdf_.pdfpt(0,xa,ev.central().muf,a.type)*
pdf_.pdfpt(1,xb,ev.central().muf,b.type);
known_pdf.emplace(ev.central().muf, result.central);
result.variations.reserve(ev.variations().size());
for(auto const & param: ev.variations()){
const double muf = param.muf;
auto cur_pdf = known_pdf.find(muf);
if(cur_pdf == known_pdf.end()){
cur_pdf = known_pdf.emplace(
muf,
pdf_.pdfpt(0,xa,muf,a.type)*pdf_.pdfpt(1,xb,muf,b.type)
).first;
}
result.variations.emplace_back(cur_pdf->second);
}
assert(result.variations.size() == ev.variations().size());
return result;
}
EventReweighter::EventFactors
EventReweighter::matrix_elements(Event const & ev) const{
assert(treat_.count(last_event_class_) > 0);
if(treat_.find(last_event_class_)->second == EventTreatment::keep){
return fixed_order_scale_ME(ev);
}
// precompute overall kinematic factor
const double ME_kin = MEt2_.tree_kin(ev.incoming(), ev.outgoing(), true);
EventFactors result;
std::unordered_map<double, double> known_ME;
result.central = MEt2_(
- pdf_.Halphas(ev.central().mur), ev.central().mur,
+ ev.central().mur,
ev.incoming(), ev.outgoing(),
true
);
known_ME.emplace(ev.central().mur, result.central);
result.variations.reserve(ev.variations().size());
for(auto const & param: ev.variations()){
const double mur = param.mur;
auto cur_ME = known_ME.find(mur);
if(cur_ME == known_ME.end()){
- const double alpha_s = pdf_.Halphas(mur);
const double ME = MEt2_.tree_param(
- alpha_s, mur, ev.incoming(), ev.outgoing()
+ mur, ev.incoming(), ev.outgoing()
)*ME_kin*MEt2_.virtual_corrections(
- alpha_s, mur, ev.incoming(), ev.outgoing()
+ mur, ev.incoming(), ev.outgoing()
);
cur_ME = known_ME.emplace(mur, ME).first;
}
result.variations.emplace_back(cur_ME->second);
}
assert(result.variations.size() == ev.variations().size());
return result;
}
double EventReweighter::tree_matrix_element(Event const & ev) const{
assert(ev.variations().empty());
assert(treat_.count(last_event_class_) > 0);
if(treat_.find(last_event_class_)->second == EventTreatment::keep){
return fixed_order_scale_ME(ev).central;
}
return MEt2_.tree(
- pdf_.Halphas(ev.central().mur), ev.central().mur,
+ ev.central().mur,
ev.incoming(), ev.outgoing(),
false
);
}
EventReweighter::EventFactors
EventReweighter::fixed_order_scale_ME(Event const & ev) const{
const int alpha_s_power = std::count_if(
begin(ev.outgoing()), end(ev.outgoing()),
[](Sparticle const & p){ return is_parton(p); }
);
EventFactors result;
result.central = pow(pdf_.Halphas(ev.central().mur), alpha_s_power);
for(auto const & var: ev.variations()){
result.variations.emplace_back(
pow(pdf_.Halphas(var.mur), alpha_s_power)
);
}
return result;
}
}
diff --git a/src/MatrixElement.cc b/src/MatrixElement.cc
index dc76c2b..58731ad 100644
--- a/src/MatrixElement.cc
+++ b/src/MatrixElement.cc
@@ -1,758 +1,763 @@
#include "RHEJ/MatrixElement.hh"
#include <CLHEP/Random/Randomize.h>
#include <CLHEP/Random/RanluxEngine.h>
#include "RHEJ/currents.hh"
#include "RHEJ/PDG_codes.hh"
#include "RHEJ/uno.hh"
#include "RHEJ/debug.hh"
namespace{
constexpr int N_C = 3;
constexpr int C_A = N_C;
constexpr double C_F = (N_C*N_C - 1.)/(2.*N_C);
constexpr double t_f = 1./2.;
constexpr double clambda = 0.2;
constexpr int n_f = 5;
constexpr double beta0 = 11./3.*C_A - 4./3.*t_f*n_f;
}
namespace RHEJ{
//cf. last line of eq. (22) in \ref Andersen:2011hs
double MatrixElement::omega0(
double alpha_s, double mur,
fastjet::PseudoJet const & q_j, double lambda
) const {
const double result = - alpha_s*N_C/M_PI*log(q_j.perp2()/(lambda*lambda));
if(! log_corr_) return result;
// use alpha_s(sqrt(q_j*lambda)), evolved to mur
return (
1. + alpha_s/(4.*M_PI)*beta0*log(mur*mur/(q_j.perp()*lambda))
)*result;
}
double MatrixElement::virtual_corrections(
- double alpha_s, double mur,
+ double mur,
std::array<Sparticle, 2> const & in,
std::vector<Sparticle> const & out
) const{
fastjet::PseudoJet const & pa = in.front().p;
#ifndef NDEBUG
fastjet::PseudoJet const & pb = in.back().p;
#endif
assert(std::is_sorted(out.begin(), out.end(), rapidity_less{}));
assert(out.size() >= 2);
assert(pa.pz() < pb.pz());
fastjet::PseudoJet q = pa - out[0].p;
size_t first_idx = 0;
size_t last_idx = out.size() - 1;
// if there is a Higgs or unordered gluon outside the extremal partons
// then it is not part of the FKL ladder and does not contribute
// to the virtual corrections
if(out.front().type == pid::Higgs || has_unob_gluon(in, out)){
q -= out[1].p;
++first_idx;
}
if(out.back().type == pid::Higgs || has_unof_gluon(in, out)){
--last_idx;
}
double exponent = 0;
+ const double alpha_s = alpha_s_(mur);
for(size_t j = first_idx; j < last_idx; ++j){
exponent += omega0(alpha_s, mur, q, clambda)*(
out[j+1].rapidity() - out[j].rapidity()
);
q -= out[j+1].p;
}
assert(
nearby(q, -1*pb)
|| out.back().type == pid::Higgs
|| has_unof_gluon(in, out)
);
return exp(exponent);
}
}
namespace {
//! Lipatov vertex for partons emitted into extremal jets
double C2Lipatov(CLHEP::HepLorentzVector qav, CLHEP::HepLorentzVector qbv, CLHEP::HepLorentzVector p1, CLHEP::HepLorentzVector p2)
{
CLHEP::HepLorentzVector temptrans=-(qav+qbv);
CLHEP::HepLorentzVector p5=qav-qbv;
CLHEP::HepLorentzVector CL=temptrans+(qav.m2()/p5.dot(p1)+2.*p5.dot(p2)/p1.dot(p2))*p1-p2*(qbv.m2()/p5.dot(p2)+2.*p5.dot(p1)/p1.dot(p2));
#if printoutput
cout << "#Fadin qa : "<<qav<<endl;
cout << "#Fadin qb : "<<qbv<<endl;
cout << "#Fadin p1 : "<<p1<<endl;
cout << "#Fadin p2 : "<<p2<<endl;
cout << "#Fadin p5 : "<<p5<<endl;
cout << "#Fadin Gauge Check : "<< CL.dot(p5)<<endl;
cout << "#Fadin C2L : "<< -CL.dot(CL)<<" "<<-CL.dot(CL)/(qav.m2()*qbv.m2())/(4./p5.perp2())<<endl;
#endif
#if 0
if (-CL.dot(CL)<0.)
// if (fabs(CL.dot(p5))>fabs(CL.dot(CL))) // not sufficient!
return 0.;
else
#endif
return -CL.dot(CL);
}
//! Lipatov vertex with soft subtraction for partons emitted into extremal jets
double C2Lipatovots(CLHEP::HepLorentzVector qav, CLHEP::HepLorentzVector qbv, CLHEP::HepLorentzVector p1, CLHEP::HepLorentzVector p2)
{
double kperp=(qav-qbv).perp();
if (kperp>clambda)
return C2Lipatov(qav, qbv, p1, p2)/(qav.m2()*qbv.m2());
else {
double Cls=(C2Lipatov(qav, qbv, p1, p2)/(qav.m2()*qbv.m2()));
double temp=Cls-4./(kperp*kperp);
// std::cout <<kperp <<" "<<temp<<" "<<4./(kperp*kperp)<<" "<<(C2Lipatov(qav, qbv, pa, pb, p1, p2)/(qav.m2()*qbv.m2()))<<std::endl;
return temp;
}
}
//! Lipatov vertex
double C2Lipatov(CLHEP::HepLorentzVector qav, CLHEP::HepLorentzVector qbv, CLHEP::HepLorentzVector pim, CLHEP::HepLorentzVector pip,CLHEP::HepLorentzVector pom, CLHEP::HepLorentzVector pop) // B
{
CLHEP::HepLorentzVector temptrans=-(qav+qbv);
CLHEP::HepLorentzVector p5=qav-qbv;
CLHEP::HepLorentzVector CL=temptrans+qav.m2()*(1./p5.dot(pip)*pip+1./p5.dot(pop)*pop)/2.-qbv.m2()*(1./p5.dot(pim)*pim+1./p5.dot(pom)*pom)/2.+(pip*(p5.dot(pim)/pip.dot(pim)+p5.dot(pom)/pip.dot(pom)) + pop*(p5.dot(pim)/pop.dot(pim)+p5.dot(pom)/pop.dot(pom)) - pim*(p5.dot(pip)/pip.dot(pim) + p5.dot(pop)/pop.dot(pim)) - pom*(p5.dot(pip)/pip.dot(pom) + p5.dot(pop)/pop.dot(pom)) )/2.;
return -CL.dot(CL);
}
//! Lipatov vertex with soft subtraction
double C2Lipatovots(CLHEP::HepLorentzVector qav, CLHEP::HepLorentzVector qbv, CLHEP::HepLorentzVector pa, CLHEP::HepLorentzVector pb, CLHEP::HepLorentzVector p1, CLHEP::HepLorentzVector p2)
{
double kperp=(qav-qbv).perp();
if (kperp>clambda)
return C2Lipatov(qav, qbv, pa, pb, p1, p2)/(qav.m2()*qbv.m2());
else {
double Cls=(C2Lipatov(qav, qbv, pa, pb, p1, p2)/(qav.m2()*qbv.m2()));
double temp=Cls-4./(kperp*kperp);
return temp;
}
}
/** Matrix element squared for tree-level current-current scattering
* @param aptype Particle a PDG ID
* @param bptype Particle b PDG ID
* @param pn Particle n Momentum
* @param pb Particle b Momentum
* @param p1 Particle 1 Momentum
* @param pa Particle a Momentum
* @returns ME Squared for Tree-Level Current-Current Scattering
*/
double ME_current(
int aptype, int bptype,
CLHEP::HepLorentzVector const & pn,
CLHEP::HepLorentzVector const & pb,
CLHEP::HepLorentzVector const & p1,
CLHEP::HepLorentzVector const & pa
){
if (aptype==21&&bptype==21) {
return jM2gg(pn,pb,p1,pa);
} else if (aptype==21&&bptype!=21) {
if (bptype > 0)
return jM2qg(pn,pb,p1,pa);
else
return jM2qbarg(pn,pb,p1,pa);
}
else if (bptype==21&&aptype!=21) { // ----- || -----
if (aptype > 0)
return jM2qg(p1,pa,pn,pb);
else
return jM2qbarg(p1,pa,pn,pb);
}
else { // they are both quark
if (bptype>0) {
if (aptype>0)
return jM2qQ(pn,pb,p1,pa);
else
return jM2qQbar(pn,pb,p1,pa);
}
else {
if (aptype>0)
return jM2qQbar(p1,pa,pn,pb);
else
return jM2qbarQbar(pn,pb,p1,pa);
}
}
throw std::logic_error("unknown particle types");
}
/** \brief Matrix element squared for tree-level current-current scattering with Higgs
* @param aptype Particle a PDG ID
* @param bptype Particle b PDG ID
* @param pn Particle n Momentum
* @param pb Particle b Momentum
* @param p1 Particle 1 Momentum
* @param pa Particle a Momentum
* @param qH t-channel momentum before Higgs
* @param qHp1 t-channel momentum after Higgs
* @returns ME Squared for Tree-Level Current-Current Scattering with Higgs
*/
double ME_Higgs_current(
int aptype, int bptype,
CLHEP::HepLorentzVector const & pn,
CLHEP::HepLorentzVector const & pb,
CLHEP::HepLorentzVector const & p1,
CLHEP::HepLorentzVector const & pa,
CLHEP::HepLorentzVector const & qH, // t-channel momentum before Higgs
CLHEP::HepLorentzVector const & qHp1, // t-channel momentum after Higgs
double mt, bool include_bottom, double mb
){
if (aptype==21&&bptype==21) // gg initial state
return MH2gg(pn,pb,p1,pa,-qHp1,-qH,mt,include_bottom,mb);
else if (aptype==21&&bptype!=21) {
if (bptype > 0)
return MH2qg(pn,pb,p1,pa,-qHp1,-qH,mt,include_bottom,mb)*4./9.;
else
return MH2qbarg(pn,pb,p1,pa,-qHp1,-qH,mt,include_bottom,mb)*4./9.;
}
else if (bptype==21&&aptype!=21) {
if (aptype > 0)
return MH2qg(p1,pa,pn,pb,-qH,-qHp1,mt,include_bottom,mb)*4./9.;
else
return MH2qbarg(p1,pa,pn,pb,-qH,-qHp1,mt,include_bottom,mb)*4./9.;
}
else { // they are both quark
if (bptype>0) {
if (aptype>0)
return MH2qQ(pn,pb,p1,pa,-qHp1,-qH,mt,include_bottom,mb)*4.*4./(9.*9.);
else
return MH2qQbar(pn,pb,p1,pa,-qHp1,-qH,mt,include_bottom,mb)*4.*4./(9.*9.);
}
else {
if (aptype>0)
return MH2qQbar(p1,pa,pn,pb,-qH,-qHp1,mt,include_bottom,mb)*4.*4./(9.*9.);
else
return MH2qbarQbar(pn,pb,p1,pa,-qHp1,-qH,mt,include_bottom,mb)*4.*4./(9.*9.);
}
}
throw std::logic_error("unknown particle types");
}
/** \brief Current matrix element squared with Higgs and unordered forward emission
* @param aptype Particle A PDG ID
* @param bptype Particle B PDG ID
* @param punof Unordered Particle Momentum
* @param pn Particle n Momentum
* @param pb Particle b Momentum
* @param p1 Particle 1 Momentum
* @param pa Particle a Momentum
* @param qH t-channel momentum before Higgs
* @param qHp1 t-channel momentum after Higgs
* @returns ME Squared with Higgs and unordered forward emission
*/
double ME_Higgs_current_unof(
int aptype, int bptype,
CLHEP::HepLorentzVector const & punof,
CLHEP::HepLorentzVector const & pn,
CLHEP::HepLorentzVector const & pb,
CLHEP::HepLorentzVector const & p1,
CLHEP::HepLorentzVector const & pa,
CLHEP::HepLorentzVector const & qH, // t-channel momentum before Higgs
CLHEP::HepLorentzVector const & qHp1, // t-channel momentum after Higgs
double mt, bool include_bottom, double mb
){
if (aptype==21&&bptype!=21) {
if (bptype > 0)
return jM2unogqHg(punof,pn,pb,p1,pa,-qHp1,-qH,mt,include_bottom,mb);
else
return jM2unogqbarHg(punof,pn,pb,p1,pa,-qHp1,-qH,mt,include_bottom,mb);
}
else { // they are both quark
if (bptype>0) {
if (aptype>0)
return jM2unogqHQ(punof,pn,pb,p1,pa,-qHp1,-qH,mt,include_bottom,mb);
else
return jM2unogqHQbar(punof,pn,pb,p1,pa,-qHp1,-qH,mt,include_bottom,mb);
}
else {
if (aptype>0)
return jM2unogqbarHQ(punof,pn,pb,p1,pa,-qHp1,-qH,mt,include_bottom,mb);
else
return jM2unogqbarHQbar(punof,pn,pb,p1,pa,-qHp1,-qH,mt,include_bottom,mb);
}
}
throw std::logic_error("unknown particle types");
}
/** \brief Current matrix element squared with Higgs and unordered backward emission
* @param aptype Particle A PDG ID
* @param bptype Particle B PDG ID
* @param pn Particle n Momentum
* @param pb Particle b Momentum
* @param punob Unordered back Particle Momentum
* @param p1 Particle 1 Momentum
* @param pa Particle a Momentum
* @param qH t-channel momentum before Higgs
* @param qHp1 t-channel momentum after Higgs
* @returns ME Squared with Higgs and unordered backward emission
*/
double ME_Higgs_current_unob(
int aptype, int bptype,
CLHEP::HepLorentzVector const & pn,
CLHEP::HepLorentzVector const & pb,
CLHEP::HepLorentzVector const & punob,
CLHEP::HepLorentzVector const & p1,
CLHEP::HepLorentzVector const & pa,
CLHEP::HepLorentzVector const & qH, // t-channel momentum before Higgs
CLHEP::HepLorentzVector const & qHp1, // t-channel momentum after Higgs
double mt, bool include_bottom, double mb
){
if (bptype==21&&aptype!=21) {
if (aptype > 0)
return jM2unobgHQg(pn,pb,punob,p1,pa,-qHp1,-qH,mt,include_bottom,mb);
else
return jM2unobgHQbarg(pn,pb,punob,p1,pa,-qHp1,-qH,mt,include_bottom,mb);
}
else { // they are both quark
if (aptype>0) {
if (bptype>0)
return jM2unobqHQg(pn,pb,punob,p1,pa,-qHp1,-qH,mt,include_bottom,mb);
else
return jM2unobqbarHQg(pn,pb,punob,p1,pa,-qHp1,-qH,mt,include_bottom,mb);
}
else {
if (bptype>0)
return jM2unobqHQbarg(pn,pb,punob,p1,pa,-qHp1,-qH,mt,include_bottom,mb);
else
return jM2unobqbarHQbarg(pn,pb,punob,p1,pa,-qHp1,-qH,mt,include_bottom,mb);
}
}
throw std::logic_error("unknown particle types");
}
CLHEP::HepLorentzVector to_HepLorentzVector(RHEJ::Sparticle const & particle){
return {particle.p.px(), particle.p.py(), particle.p.pz(), particle.p.E()};
}
} // RHEJ namespace
namespace RHEJ{
MatrixElement::MatrixElement(
+ std::function<double (double)> alpha_s,
fastjet::JetDefinition jet_def, double jetptmin,
bool log_corr,
HiggsCouplingSettings Hgg_settings
):
+ alpha_s_{std::move(alpha_s)},
Hgg_settings_{Hgg_settings},
jet_def_{jet_def},
jetptmin_{jetptmin},
log_corr_{log_corr}
{}
double MatrixElement::operator()(
- double alpha_s, double mur,
+ double mur,
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> const & outgoing,
bool check_momenta
) const {
return tree(
- alpha_s, mur,
+ mur,
incoming, outgoing,
check_momenta
)*virtual_corrections(
- alpha_s, mur,
+ mur,
incoming, outgoing
);
}
double MatrixElement::tree_kin(
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> const & outgoing,
bool check_momenta
) const {
assert(
std::is_sorted(
incoming.begin(), incoming.end(),
[](Sparticle o1, Sparticle o2){return o1.p.pz()<o2.p.pz();}
)
);
assert(std::is_sorted(outgoing.begin(), outgoing.end(), rapidity_less{}));
auto AWZH_boson = std::find_if(
begin(outgoing), end(outgoing),
[](Sparticle const & p){return is_AWZH_boson(p);}
);
if(AWZH_boson == end(outgoing)){
return tree_kin_jets(incoming, outgoing, check_momenta);
}
switch(AWZH_boson->type){
case pid::Higgs: {
- constexpr double alpha_s_mH = 0.113559;
+ static constexpr double mH = 125.;
+ const double alpha_s_mH = alpha_s_(mH);
return alpha_s_mH*alpha_s_mH/(256.*pow(M_PI, 5))*tree_kin_Higgs(
incoming, outgoing, check_momenta
);
}
// TODO
case pid::Wp:
case pid::Wm:
case pid::photon:
case pid::Z:
default:
throw std::logic_error("Emission of boson of unsupported type");
}
}
namespace{
constexpr int extremal_jet_idx = 1;
constexpr int no_extremal_jet_idx = 0;
bool treat_as_extremal(Sparticle const & parton){
return parton.p.user_index() == extremal_jet_idx;
}
template<class InputIterator>
double FKL_ladder_weight(
InputIterator begin_gluon, InputIterator end_gluon,
CLHEP::HepLorentzVector const & q0,
CLHEP::HepLorentzVector const & pa, CLHEP::HepLorentzVector const & pb,
CLHEP::HepLorentzVector const & p1, CLHEP::HepLorentzVector const & pn
){
double wt = 1;
auto qi = q0;
for(auto gluon_it = begin_gluon; gluon_it != end_gluon; ++gluon_it){
assert(gluon_it->type == pid::gluon);
const auto g = to_HepLorentzVector(*gluon_it);
const auto qip1 = qi - g;
if(treat_as_extremal(*gluon_it)){
wt *= C2Lipatovots(qip1, qi, pa, pb)*C_A;
} else{
wt *= C2Lipatovots(qip1, qi, pa, pb, p1, pn)*C_A;
}
qi = qip1;
}
return wt;
}
} // anonymous namespace
std::vector<Sparticle> MatrixElement::tag_extremal_jet_partons(
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> out_partons, bool check_momenta
) const{
if(!check_momenta){
for(auto & parton: out_partons){
parton.p.set_user_index(no_extremal_jet_idx);
}
return out_partons;
}
fastjet::ClusterSequence cs(to_PseudoJet(out_partons), jet_def_);
const auto jets = sorted_by_rapidity(cs.inclusive_jets(jetptmin_));
assert(jets.size() >= 2);
auto most_backward = begin(jets);
auto most_forward = end(jets) - 1;
// skip jets caused by unordered emission
if(has_unob_gluon(incoming, out_partons)){
assert(jets.size() >= 3);
++most_backward;
}
else if(has_unof_gluon(incoming, out_partons)){
assert(jets.size() >= 3);
--most_forward;
}
const auto extremal_jet_indices = cs.particle_jet_indices(
{*most_backward, *most_forward}
);
assert(extremal_jet_indices.size() == out_partons.size());
for(size_t i = 0; i < out_partons.size(); ++i){
assert(RHEJ::is_parton(out_partons[i]));
const int idx = (extremal_jet_indices[i]>=0)?
extremal_jet_idx:
no_extremal_jet_idx;
out_partons[i].p.set_user_index(idx);
}
return out_partons;
}
double MatrixElement::tree_kin_jets(
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> partons,
bool check_momenta
) const {
partons = tag_extremal_jet_partons(incoming, partons, check_momenta);
if(has_unob_gluon(incoming, partons) || has_unof_gluon(incoming, partons)){
throw std::logic_error("unordered emission not implemented for pure jets");
}
const auto pa = to_HepLorentzVector(incoming[0]);
const auto pb = to_HepLorentzVector(incoming[1]);
const auto p1 = to_HepLorentzVector(partons.front());
const auto pn = to_HepLorentzVector(partons.back());
return ME_current(
incoming[0].type, incoming[1].type,
pn, pb, p1, pa
)/(4*(N_C*N_C - 1))*FKL_ladder_weight(
begin(partons) + 1, end(partons) - 1,
pa - p1, pa, pb, p1, pn
);
}
double MatrixElement::tree_kin_Higgs(
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> const & outgoing,
bool check_momenta
) const {
if(has_uno_gluon(incoming, outgoing)){
return tree_kin_Higgs_between(incoming, outgoing, check_momenta);
}
if(outgoing.front().type == pid::Higgs){
return tree_kin_Higgs_first(incoming, outgoing, check_momenta);
}
if(outgoing.back().type == pid::Higgs){
return tree_kin_Higgs_last(incoming, outgoing, check_momenta);
}
return tree_kin_Higgs_between(incoming, outgoing, check_momenta);
}
double MatrixElement::MH2_forwardH(
RHEJ::ParticleID id,
CLHEP::HepLorentzVector p1out, CLHEP::HepLorentzVector p1in,
CLHEP::HepLorentzVector p2out, CLHEP::HepLorentzVector p2in,
CLHEP::HepLorentzVector pH,
double t1, double t2
) const{
ignore(p2out, p2in);
const double shat = p1in.invariantMass2(p2in);
assert(RHEJ::is_parton(id));
if(id != RHEJ::pid::gluon){
return 9./2.*shat*shat*C2qHqm(p1in,p1out,pH)/(t1*t2);
}
// gluon case
#ifdef RHEJ_BUILD_WITH_QCDLOOP
if(!Hgg_settings_.use_impact_factors){
return C_A/C_F*1./(16*M_PI*M_PI)*t1/t2*MH2gq_outsideH(
p1out, p1in, p2out, p2in, pH,
Hgg_settings_.mt, Hgg_settings_.include_bottom,
Hgg_settings_.mb
);
}
#endif
return 9./2.*shat*shat*(
C2gHgp(p1in,p1out,pH) + C2gHgm(p1in,p1out,pH)
)/(t1*t2);
}
double MatrixElement::tree_kin_Higgs_first(
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> const & outgoing,
bool check_momenta
) const {
assert(outgoing.front().type == pid::Higgs);
const auto pH = to_HepLorentzVector(outgoing.front());
const auto partons = tag_extremal_jet_partons(
incoming,
std::vector<Sparticle>(begin(outgoing) + 1, end(outgoing)),
check_momenta
);
const auto pa = to_HepLorentzVector(incoming[0]);
const auto pb = to_HepLorentzVector(incoming[1]);
const auto p1 = to_HepLorentzVector(partons.front());
const auto pn = to_HepLorentzVector(partons.back());
const auto q0 = pa - p1 - pH;
const double t1 = q0.m2();
const double t2 = (pn - pb).m2();
double wt = MH2_forwardH(
incoming[0].type, p1, pa, pn, pb, pH,
t1, t2
)*FKL_ladder_weight(
begin(partons) + 1, end(partons) - 1,
q0, pa, pb, p1, pn
);
for(auto const & inc: incoming){
if(inc.type != pid::gluon) wt *= C_F/C_A;
}
return wt;
}
double MatrixElement::tree_kin_Higgs_last(
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> const & outgoing,
bool check_momenta
) const {
assert(outgoing.back().type == pid::Higgs);
const auto pH = to_HepLorentzVector(outgoing.back());
const auto partons = tag_extremal_jet_partons(
incoming,
std::vector<Sparticle>(begin(outgoing), end(outgoing) - 1),
check_momenta
);
const auto pa = to_HepLorentzVector(incoming[0]);
const auto pb = to_HepLorentzVector(incoming[1]);
auto p1 = to_HepLorentzVector(partons.front());
const auto pn = to_HepLorentzVector(partons.back());
auto q0 = pa - p1;
const double t1 = q0.m2();
const double t2 = (pn + pH - pb).m2();
double wt = MH2_forwardH(
incoming[1].type, pn, pb, p1, pa, pH,
t2, t1
)*FKL_ladder_weight(
begin(partons) + 1, end(partons) - 1,
q0, pa, pb, p1, pn
);
for(auto const & inc: incoming){
if(inc.type != pid::gluon) wt *= C_F/C_A;
}
return wt;
}
double MatrixElement::tree_kin_Higgs_between(
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> const & outgoing,
bool check_momenta
) const {
const auto the_Higgs = std::find_if(
begin(outgoing), end(outgoing),
[](Sparticle const & s){ return s.type == pid::Higgs; }
);
assert(the_Higgs != end(outgoing));
const auto pH = to_HepLorentzVector(*the_Higgs);
std::vector<Sparticle> partons(begin(outgoing), the_Higgs);
partons.insert(end(partons), the_Higgs + 1, end(outgoing));
partons = tag_extremal_jet_partons(incoming, partons, check_momenta);
const auto pa = to_HepLorentzVector(incoming[0]);
const auto pb = to_HepLorentzVector(incoming[1]);
auto p1 = to_HepLorentzVector(
partons[has_unob_gluon(incoming, outgoing)?1:0]
);
auto pn = to_HepLorentzVector(
partons[partons.size() - (has_unof_gluon(incoming, outgoing)?2:1)]
);
auto first_after_Higgs = begin(partons) + (the_Higgs-begin(outgoing));
assert(
(first_after_Higgs == end(partons) && has_unob_gluon(incoming, outgoing))
|| first_after_Higgs->rapidity() >= the_Higgs->rapidity()
);
assert(
(first_after_Higgs == begin(partons) && has_unof_gluon(incoming, outgoing))
|| (first_after_Higgs-1)->rapidity() <= the_Higgs->rapidity()
);
// always treat the Higgs as if it were in between the extremal FKL partons
if(first_after_Higgs == begin(partons)) ++first_after_Higgs;
else if(first_after_Higgs == end(partons)) --first_after_Higgs;
// t-channel momentum before Higgs
auto qH = pa;
for(auto parton_it = begin(partons); parton_it != first_after_Higgs; ++parton_it){
qH -= to_HepLorentzVector(*parton_it);
}
auto q0 = pa - p1;
auto begin_ladder = begin(partons) + 1;
auto end_ladder = end(partons) - 1;
double current_factor;
if(has_unob_gluon(incoming, outgoing)){
current_factor = 9./2.*ME_Higgs_current_unob(
incoming[0].type, incoming[1].type,
pn, pb, to_HepLorentzVector(partons.front()), p1, pa, qH, qH - pH,
Hgg_settings_.mt, Hgg_settings_.include_bottom, Hgg_settings_.mb
);
const auto p_unob = to_HepLorentzVector(partons.front());
q0 -= p_unob;
p1 += p_unob;
++begin_ladder;
}
else if(has_unof_gluon(incoming, outgoing)){
current_factor = 9./2.*ME_Higgs_current_unof(
incoming[0].type, incoming[1].type,
to_HepLorentzVector(partons.back()), pn, pb, p1, pa, qH, qH - pH,
Hgg_settings_.mt, Hgg_settings_.include_bottom, Hgg_settings_.mb
);
pn += to_HepLorentzVector(partons.back());
--end_ladder;
}
else{
current_factor = ME_Higgs_current(
incoming[0].type, incoming[1].type,
pn, pb, p1, pa, qH, qH - pH,
Hgg_settings_.mt, Hgg_settings_.include_bottom, Hgg_settings_.mb
);
}
const double ladder_factor = FKL_ladder_weight(
begin_ladder, first_after_Higgs,
q0, pa, pb, p1, pn
)*FKL_ladder_weight(
first_after_Higgs, end_ladder,
qH - pH, pa, pb, p1, pn
);
return current_factor*9./8.*ladder_factor;
}
double MatrixElement::tree_param_partons(
double alpha_s, double mur,
std::vector<Sparticle> const & partons
) const{
const double gs2 = 4.*M_PI*alpha_s;
double wt = std::pow(gs2, partons.size());
if(log_corr_){
// use alpha_s(q_perp), evolved to mur
assert(partons.size() >= 2);
for(size_t i = 1; i < partons.size()-1; ++i){
wt *= 1 + alpha_s/(2*M_PI)*beta0*log(mur/partons[i].p.perp());
}
}
return wt;
}
double MatrixElement::tree_param(
- double alpha_s, double mur,
+ double mur,
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> const & outgoing
) const{
+ const double alpha_s = alpha_s_(mur);
if(has_unob_gluon(incoming, outgoing)){
return 4*M_PI*alpha_s*tree_param_partons(
alpha_s, mur, filter_partons({begin(outgoing) + 1, end(outgoing)})
);
}
if(has_unof_gluon(incoming, outgoing)){
return 4*M_PI*alpha_s*tree_param_partons(
alpha_s, mur, filter_partons({begin(outgoing), end(outgoing) - 1})
);
}
return tree_param_partons(alpha_s, mur, filter_partons(outgoing));
}
double MatrixElement::tree(
- double alpha_s, double mur,
+ double mur,
std::array<Sparticle, 2> const & incoming,
std::vector<Sparticle> const & outgoing,
bool check_momenta
) const {
- return tree_param(alpha_s, mur, incoming, outgoing)*tree_kin(
+ return tree_param(mur, incoming, outgoing)*tree_kin(
incoming, outgoing, check_momenta
);
}
} // namespace RHEJ
diff --git a/src/config.cc b/src/config.cc
index e7e9f5e..2468110 100644
--- a/src/config.cc
+++ b/src/config.cc
@@ -1,606 +1,458 @@
#include "RHEJ/config.hh"
#include "RHEJ/exceptions.hh"
#include <set>
#include <string>
#include <vector>
#include <iostream>
#include <stdexcept>
namespace RHEJ{
- 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 == "hepmc3") return FileFormat::HepMC;
- throw std::invalid_argument{
- "Can't determine format for output file " + filename
- };
- }
-
-}
-
-namespace YAML {
-
- std::string to_string(NodeType::value type){
- switch (type) {
- case NodeType::Null: return "Null";
- case NodeType::Scalar: return "Scalar";
- case NodeType::Sequence: return "Sequence";
- case NodeType::Map: return "Map";
- case NodeType::Undefined: return "Undefined";
- default:
- throw std::logic_error{"Missing case in switch statement"};
- }
- }
-
- template<>
- struct convert<RHEJ::OutputFile> {
- static Node encode(RHEJ::OutputFile const & outfile) {
- Node node;
- node[to_string(outfile.format)] = outfile.name;
- return node;
- }
- static bool decode(Node const & node, RHEJ::OutputFile & out) {
- switch(node.Type()){
- case NodeType::Map: {
- YAML::const_iterator it = node.begin();
- out.format = RHEJ::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 = RHEJ::format_from_suffix(out.name);
- return true;
- default:
- return false;
- }
- }
- };
-}
-
-namespace RHEJ{
namespace{
- static const std::set<std::string> known = {
- {"trials"},
- {"min extparton pt"}, {"max ext soft pt fraction"},
- {"resummation jets"}, {"fixed order jets"},
- {"FKL"}, {"non-FKL"}, {"unordered"},
- {"log correction"},
- {"event output"},
- {"analysis"},
- {"unweight"}, {"RanLux init"},
- {"scales"}, {"scale factors"}, {"max scale ratio"},
- {"Higgs coupling"}
- };
-
- static const std::set<std::string> known_jet = {
- {"min pt"}, {"algorithm"}, {"R"}
- };
-
- template<typename T>
- std::string type_string();
-
- template<>
- std::string type_string<std::vector<double>>(){
- return "array of doubles";
- }
-
- template<>
- std::string type_string<std::vector<std::string>>(){
- return "array of strings";
- }
-
- template<>
- std::string type_string<int>(){
- return "integer";
- }
-
- template<>
- std::string type_string<double>(){
- return "double";
- }
-
- template<>
- std::string type_string<bool>(){
- return "bool";
- }
-
- template<>
- std::string type_string<std::string>(){
- return "string";
- }
-
- template<>
- std::string type_string<std::vector<OutputFile>>(){
- return "array of output files";
+ //! 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-FKL",
+ "scales", "scale factors", "max scale ratio",
+ "log correction", "unweight", "event output", "analysis",
+ "RanLux init"
+ };
+ // 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] = "";
+ }
+ 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;
}
- // helper struct to allow partial template specialisation
- template<typename T>
- struct existing_as_helper{
-
- static T as(YAML::Node const & config, std::string const & entry){
- assert(config[entry]);
- try{
- return config[entry].as<T>();
- }
- catch(std::exception const &){
- throw invalid_type{
- "Entry " + entry + " is not of type " + type_string<T>()
- };
- }
- }
+ }
+ ParticleID to_ParticleID(std::string const & name){
+ using namespace RHEJ::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;
+ }
- template<>
- struct existing_as_helper<fastjet::JetAlgorithm>{
+ namespace detail{
+ void set_from_yaml(fastjet::JetAlgorithm & setting, YAML::Node const & yaml){
+ setting = to_JetAlgorithm(yaml.as<std::string>());
+ }
- static fastjet::JetAlgorithm as(
- YAML::Node const & config, std::string const & entry
- ){
- assert(config[entry]);
- const std::string algo_name =
- existing_as_helper<std::string>::as(config, entry);
- return to_JetAlgorithm(algo_name);
- }
+ 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>());
+ }
+ }
- template<>
- struct existing_as_helper<EventTreatment>{
+ 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;
+ }
- static EventTreatment as(
- YAML::Node const & config, std::string const & entry
- ){
- assert(config[entry]);
- const std::string name =
- existing_as_helper<std::string>::as(config, entry);
- return to_EventTreatment(name);
+ HiggsCouplingSettings get_Higgs_coupling(
+ YAML::Node const & node,
+ std::string const & entry
+ ){
+ assert(node);
+ static constexpr double mt_max = 2e4;
+#ifndef RHEJ_BUILD_WITH_QCDLOOP
+ if(node[entry]){
+ throw std::invalid_argument{
+ "Higgs coupling settings require building Reversed HEJ "
+ "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;
+ }
- template<>
- struct existing_as_helper<OutputFile>{
+ std::string extract_suffix(std::string const & filename){
+ size_t separator = filename.rfind('.');
+ if(separator == filename.npos) return {};
+ return filename.substr(separator + 1);
+ }
- static OutputFile as(
- YAML::Node const & config, std::string const & entry
- ){
- assert(config[entry]);
- YAML::convert<RHEJ::OutputFile> converter{};
- OutputFile out;
- if(converter.decode(config[entry], out)) return out;
- throw std::invalid_argument{
- "Bad output file setting: " + config[entry].as<std::string>()
- };
- }
+ FileFormat format_from_suffix(std::string const & filename){
+ const std::string suffix = extract_suffix(filename);
+ if(suffix == "lhe") return FileFormat::Les_Houches;
+ if(suffix == "hepmc3") return FileFormat::HepMC;
+ throw std::invalid_argument{
+ "Can't determine format for output file " + filename
};
+ }
- template<typename T>
- struct existing_as_helper<std::vector<T>>{
-
- static std::vector<T> as(
- YAML::Node const & config, std::string const & entry
- ){
- assert(config[entry]);
- // special case: treat a single value like a vector with one element
- if(config[entry].IsScalar()){
- return {existing_as_helper<T>::as(config, entry)};
- }
+ void assert_all_options_known(
+ YAML::Node const & conf, YAML::Node const & supported
+ ){
+ if(!conf.IsMap()) return;
+ 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
+ */
+ if(name != "analysis"){
try{
- return config[entry].as<std::vector<T>>();
+ assert_all_options_known(conf[name], supported[name]);
}
- catch(std::exception const & p){
- std::cerr << "Error: " << p.what() << '\n';
- throw invalid_type{
- "Entry " + entry + " is not of type " + type_string<std::vector<T>>()
- };
+ catch(unknown_option const & ex){
+ throw unknown_option{name + ": " + ex.what()};
}
}
-
- };
-
- template<typename T>
- T existing_as(YAML::Node const & config, std::string const & entry){
- return existing_as_helper<T>::as(config, entry);
}
+ }
- template<typename T>
- T as(YAML::Node const & config, std::string const & entry);
-
- template<>
- JetParameters existing_as<JetParameters>(
- YAML::Node const & config, std::string const & entry
- ){
- assert(config[entry]);
- YAML::Node const & jet_config = config[entry];
- JetParameters result;
- for(auto const & opt: jet_config){
- auto const & opt_name = opt.first.as<std::string>();
- if(! known_jet.count(opt_name)){
- throw unknown_option{"Unknown option " + opt_name};
- }
- }
- try{
- result.def = fastjet::JetDefinition{
- as<fastjet::JetAlgorithm>(jet_config, "algorithm"),
- as<double>(jet_config, "R")
- };
- result.min_pt = as<double>(jet_config, "min pt");
- return result;
- }
- catch(missing_option const & exc){
- throw missing_option{"In entry " + entry + ":\n" + exc.what()};
- }
- catch(invalid_type const & exc){
- throw invalid_type{"In entry " + entry + ":\n" + exc.what()};
- }
- }
+}
- template<typename T>
- T as(YAML::Node const & config, std::string const & entry){
- if(!config[entry]){
- throw missing_option{"No entry for option " + entry};
- }
- return existing_as<T>(config, entry);
- }
+namespace YAML {
- template<typename T>
- optional<T> as_optional(YAML::Node const & config, std::string const & entry){
- if(!config[entry]) return {};
- return {existing_as<T>(config, entry)};
- }
+ Node convert<RHEJ::OutputFile>::encode(RHEJ::OutputFile const & outfile) {
+ Node node;
+ node[to_string(outfile.format)] = outfile.name;
+ return node;
+ };
- template<>
- HiggsCouplingSettings as<HiggsCouplingSettings>(
- YAML::Node const & config, std::string const & entry
- ){
- static constexpr double mt_max = 2e4;
- HiggsCouplingSettings settings;
- if(!config[entry]) return settings;
-#ifndef RHEJ_BUILD_WITH_QCDLOOP
- throw std::invalid_argument{
- "Higgs coupling settings require building Reversed HEJ "
- "with QCDloop support"
- };
-#endif
- YAML::Node const & node = config[entry];
- if(!node.IsMap()){
- throw std::invalid_argument{
- "entry '" + entry + "' is a " + to_string(node.Type())
- + ", not a Map of properties"
- };
- }
- for(auto const & property: node){
- const auto key = property.first.as<std::string>();
- auto const & value = property.second;
- if(key == "mt") settings.mt = value.as<double>();
- else if(key == "mb") settings.mb = value.as<double>();
- else if(key == "include bottom"){
- settings.include_bottom = value.as<bool>();
- }
- else if(key == "use impact factors"){
- settings.use_impact_factors = value.as<bool>();
- }
- else{
- throw unknown_option{"Unknown option '" + key + "'"};
- }
- }
- 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;
+ bool convert<RHEJ::OutputFile>::decode(Node const & node, RHEJ::OutputFile & out) {
+ switch(node.Type()){
+ case NodeType::Map: {
+ YAML::const_iterator it = node.begin();
+ out.format = RHEJ::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 = RHEJ::format_from_suffix(out.name);
+ return true;
+ default:
+ return false;
}
+ }
+}
+namespace RHEJ{
+ namespace{
void update_fixed_order_jet_parameters(
JetParameters & fixed_order_jets, YAML::Node const & yaml
){
- YAML::Node const & yaml_jets = yaml["fixed order jets"];
- assert(yaml_jets);
- for(auto const & opt: yaml_jets){
- auto const & opt_name = opt.first.as<std::string>();
- if(! known_jet.count(opt_name)){
- throw unknown_option{"Unknown option '" + opt_name + "'"};
- }
- }
- try{
- auto algo = as_optional<fastjet::JetAlgorithm>(yaml_jets, "algorithm");
- if(algo){
- fixed_order_jets.def = fastjet::JetDefinition{
- *algo, fixed_order_jets.def.R()
- };
- }
- auto R = as_optional<double>(yaml_jets, "R");
- if(R){
- fixed_order_jets.def = fastjet::JetDefinition{
- fixed_order_jets.def.jet_algorithm(), *R
- };
- }
- auto min_pt = as_optional<double>(yaml_jets, "min pt");
- if(min_pt) fixed_order_jets.min_pt = *min_pt;
- }
- catch(missing_option const & exc){
- throw missing_option{
- std::string{"In entry fixed order jets:\n"} + exc.what()
- };
- }
- catch(invalid_type const & exc){
- throw invalid_type{
- std::string{"In entry fixed order jets:\n"} + exc.what()
- };
- }
+ 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;
}
// simple (as in non-composite) scale functions
/**
* An example for a simple scale function would be Ht,
* Ht/2 is then composite (take Ht and then divide by 2)
*/
std::unique_ptr<ScaleFunction> make_simple_ScaleFunction_ptr(
std::string const & scale_fun
){
using ret_type = std::unique_ptr<ScaleFunction>;
assert(
scale_fun.empty() ||
(!std::isspace(scale_fun.front()) && !std::isspace(scale_fun.back()))
);
if(scale_fun == "input") return ret_type{new InputScales{}};
if(scale_fun == "Ht") return ret_type{new Ht{}};
if(scale_fun == "max jet pperp") return ret_type{new MaxJetPperp{}};
if(scale_fun == "jet invariant mass"){
return ret_type{new JetInvariantMass{}};
}
try{
const double scale = to_double(scale_fun);
return ret_type{new 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;
}
std::unique_ptr<ScaleFunction> make_ScaleFunction_ptr(
std::string const & scale_fun
){
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 make_simple_ScaleFunction_ptr(scale_fun);
}
const std::string first = trim_back(std::string{scale_fun, 0, delim});
const std::string second = trim_front(std::string{scale_fun, delim+1});
double factor;
std::unique_ptr<ScaleFunction> fun;
if(scale_fun[delim] == '/'){
factor = 1/to_double(second);
fun = make_simple_ScaleFunction_ptr(first);
}
else{
assert(scale_fun[delim] == '*');
try{
factor = to_double(second);
fun = make_simple_ScaleFunction_ptr(first);
}
catch(std::invalid_argument const &){
factor = to_double(first);
fun = make_simple_ScaleFunction_ptr(second);
}
}
assert(fun != nullptr);
return std::unique_ptr<ScaleFunction>{
new Product{factor, std::move(fun)}
};
}
- ScaleConfig to_ScaleConfig(YAML::Node const & yaml){
- ScaleConfig config;
- const auto scales = as<std::vector<std::string>>(yaml, "scales");
- config.scales.reserve(scales.size());
- std::transform(
- begin(scales), end(scales), std::back_inserter(config.scales),
- make_ScaleFunction_ptr
- );
- auto scale_factors = as_optional<std::vector<double>>(
- yaml, "scale factors"
- );
- if(scale_factors) config.scale_factors = std::move(*scale_factors);
- const auto max_scale_ratio = as_optional<double>(yaml, "max scale ratio");
- static_assert(
- std::numeric_limits<double>::has_infinity, "infinity not supported"
- );
- config.max_scale_ratio = max_scale_ratio?
- *max_scale_ratio:
- std::numeric_limits<double>::infinity()
- ;
- return config;
- }
-
EventTreatMap get_event_treatment(
YAML::Node const & yaml
){
using namespace event_class;
EventTreatMap treat {
{no_2_jets, EventTreatment::discard},
- {bad_final_state, EventTreatment::discard}
+ {bad_final_state, EventTreatment::discard},
+ {FKL, EventTreatment::reweight},
+ {unob, EventTreatment::keep},
+ {unof, EventTreatment::keep},
+ {nonFKL, EventTreatment::keep}
};
- treat.emplace(FKL, as<EventTreatment>(yaml, "FKL"));
- const auto unordered_treatment = as<EventTreatment>(yaml, "unordered");
- treat.emplace(unordered_forward, unordered_treatment);
- treat.emplace(unordered_backward, unordered_treatment);
- treat.emplace(nonFKL, as<EventTreatment>(yaml, "non-FKL"));
+ 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(nonFKL), yaml, "non-FKL");
if(treat[nonFKL] == EventTreatment::reweight){
throw std::invalid_argument{"Cannot reweight non-FKL events"};
}
return treat;
}
+
Config to_Config(YAML::Node const & yaml){
- for(auto const & opt: yaml){
- auto const & opt_name = opt.first.as<std::string>();
- if(! known.count(opt_name)){
- throw unknown_option{"Unknown option '" + opt_name + "'"};
- }
+ 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 = as<JetParameters>(yaml, "resummation jets");
+ config.resummation_jets = get_jet_parameters(yaml, "resummation jets");
config.fixed_order_jets = config.resummation_jets;
- if(yaml["fixed order jets"]){
- update_fixed_order_jet_parameters(config.fixed_order_jets, yaml);
- }
- config.min_extparton_pt = as<double>(yaml, "min extparton pt");
- config.max_ext_soft_pt_fraction = yaml["max ext soft pt fraction"]?
- as<double>(yaml, "max ext soft pt fraction"):
- std::numeric_limits<double>::infinity()
- ;
- config.trials = as<int>(yaml, "trials");
- config.log_correction = as<bool>(yaml, "log correction");
- config.unweight = as<bool>(yaml, "unweight");
+ 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);
- if(yaml["event output"]){
- config.output = as<std::vector<OutputFile>>(yaml, "event output");
- }
- config.RanLux_init = as_optional<std::string>(yaml, "RanLux init");
- config.analysis_parameters = yaml["analysis"]?yaml["analysis"]:YAML::Node{};
+ set_from_yaml_if_defined(config.output, yaml, "event output");
+ set_from_yaml_if_defined(config.RanLux_init, yaml, "RanLux init");
+ set_from_yaml_if_defined(config.analysis_parameters, yaml, "analysis");
config.scale_gen = ScaleGenerator{to_ScaleConfig(yaml)};
- config.Higgs_coupling = as<HiggsCouplingSettings>(yaml, "Higgs coupling");
+ config.Higgs_coupling = get_Higgs_coupling(yaml, "Higgs coupling");
return config;
}
} // anonymous namespace
+ ScaleConfig to_ScaleConfig(YAML::Node const & yaml){
+ ScaleConfig config;
+ std::vector<std::string> scales;
+ set_from_yaml(scales, yaml, "scales");
+ config.scales.reserve(scales.size());
+ std::transform(
+ begin(scales), end(scales), std::back_inserter(config.scales),
+ make_ScaleFunction_ptr
+ );
+ set_from_yaml_if_defined(config.scale_factors, yaml, "scale factors");
+ config.max_scale_ratio = std::numeric_limits<double>::infinity();
+ set_from_yaml_if_defined(config.max_scale_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;
}
}
Event ScaleGenerator::operator()(Event ev) const{
assert(ev.variations().empty());
assert(! settings_.scales.empty());
auto const & scales = settings_.scales;
auto const & scale_factors = settings_.scale_factors;
const double max_scale_ratio = settings_.max_scale_ratio;
// check if we are doing scale variation at all
if(scales.size() == 1 && scale_factors.empty()){
ev.central() = (*scales.front())(ev);
return ev;
}
for(auto && base_scale: scales){
const EventParameters cur_base = (*base_scale)(ev);
ev.variations().emplace_back(cur_base);
//multiplicative scale variation
for(double mur_factor: scale_factors){
const double mur = cur_base.mur*mur_factor;
for(double muf_factor: scale_factors){
if(muf_factor == 1. && mur_factor == 1.) continue;
const double muf = cur_base.muf*muf_factor;
if(mur/muf < 1/max_scale_ratio || mur/muf > max_scale_ratio) continue;
ev.variations().emplace_back(
EventParameters{mur, muf, cur_base.weight}
);
}
}
};
ev.central() = (*scales.front())(ev);
return ev;
}
}
diff --git a/src/main.cc b/src/main.cc
index 5d67c5f..763cf65 100644
--- a/src/main.cc
+++ b/src/main.cc
@@ -1,144 +1,147 @@
/**
* Name: main.cc
* Authors: Tuomas Hapola, Andreas Maier <andreas.maier@durham.ac.uk>
*
*/
#include <fstream>
#include <algorithm>
#include <memory>
#include <chrono>
#include <iostream>
#include "yaml-cpp/yaml.h"
#include "LHEF/LHEF.h"
#include "RHEJ/CombinedEventWriter.hh"
#include "RHEJ/get_analysis.hh"
#include "RHEJ/utility.hh"
#include "RHEJ/EventReweighter.hh"
#include "RHEJ/config.hh"
#include "RHEJ/stream.hh"
int event_number(std::string const & record){
size_t start = record.rfind("Number of Events");
start = record.find_first_of("123456789", start);
if(start == std::string::npos) {
throw std::invalid_argument("no event number record found");
}
const size_t end = record.find_first_not_of("0123456789", start);
return std::stoi(record.substr(start, end - start));
}
RHEJ::Config load_config(char const * filename){
try{
return RHEJ::load_config(filename);
}
catch(std::exception const & exc){
std::cerr << "Error: " << exc.what() << '\n';
std::exit(EXIT_FAILURE);
}
}
std::unique_ptr<RHEJ::Analysis> get_analysis(
YAML::Node const & parameters
){
try{
return RHEJ::get_analysis(parameters);
}
catch(std::exception const & exc){
std::cerr << "Failed to load analysis: " << exc.what() << '\n';
std::exit(EXIT_FAILURE);
}
}
int main(int argn, char** argv) {
using clock = std::chrono::system_clock;
if (argn < 3) {
std::cerr << "\n# Usage:\n./HEJ config_file input_file\n\n";
return EXIT_FAILURE;
}
const auto start_time = clock::now();
// read configuration
const RHEJ::Config config = load_config(argv[1]);
RHEJ::istream in{argv[2]};
std::unique_ptr<RHEJ::Analysis> analysis = get_analysis(
config.analysis_parameters
);
assert(analysis != nullptr);
LHEF::Reader reader{in};
RHEJ::CombinedEventWriter writer{config.output, reader.heprup};
double global_reweight = 1.;
int max_events = std::numeric_limits<int>::max();
if(argn > 3){
max_events = std::stoi(argv[3]);
const int input_events = event_number(reader.headerBlock);
global_reweight = input_events/static_cast<double>(max_events);
std::cout << "processing " << max_events
<< " out of " << input_events << " events\n";
}
RHEJ::EventReweighter rhej{
reader.heprup,
config.resummation_jets.def, config.resummation_jets.min_pt,
config.treat,
config.min_extparton_pt, config.max_ext_soft_pt_fraction,
config.log_correction,
config.Higgs_coupling
};
if(config.RanLux_init){
RHEJ::PhaseSpacePoint::reset_ranlux(*config.RanLux_init);
}
int nevent = 0;
std::array<int, RHEJ::event_class::last_class + 1>
nevent_class{0}, nfailed_class{0};
// Loop over the events in the inputfile
while(reader.readEvent()){
// reweight events so that the total cross section is conserved
reader.hepeup.setWeight(0, global_reweight * reader.hepeup.weight());
if(nevent == max_events) break;
nevent++;
if (nevent % 10000 == 0)
std::cout << "event number " << nevent << std::endl;
// calculate rHEJ weight
RHEJ::Event FO_event{
RHEJ::UnclusteredEvent{reader.hepeup},
config.fixed_order_jets.def, config.fixed_order_jets.min_pt,
};
auto resummed_events = rhej.reweight(
FO_event,
config.trials, config.scale_gen
);
const auto ev_class = rhej.last_event_class();
++nevent_class[ev_class];
if(resummed_events.empty()) ++nfailed_class[ev_class];
for(auto const & ev: resummed_events){
- analysis->fill(ev);
- writer.write(ev);
+ //TODO: move pass_cuts to after phase space point generation
+ if(analysis->pass_cuts(ev)){
+ analysis->fill(ev);
+ writer.write(ev);
+ }
}
} // main event loop
using namespace RHEJ::event_class;
std::cout<< "Events processed: " << nevent << '\n';
for(size_t ev_class = first_class; ev_class <= last_class; ++ev_class){
std::cout << '\t' << names[ev_class] << ": " << nevent_class[ev_class]
<< ", failed to reconstruct " << nfailed_class[ev_class]
<< '\n';
}
std::chrono::duration<double> run_time = (clock::now() - start_time);
std::cout << "\nTask Runtime: " << run_time.count() << " seconds.\n";
}
diff --git a/t/test_ME_h_3j.cc b/t/test_ME_h_3j.cc
index bf703a3..dffaa38 100644
--- a/t/test_ME_h_3j.cc
+++ b/t/test_ME_h_3j.cc
@@ -1,87 +1,93 @@
#include "LHEF/LHEF.h"
#include "RHEJ/MatrixElement.hh"
#include "RHEJ/Event.hh"
struct Event{
std::array<RHEJ::Sparticle, 2> incoming;
std::vector<RHEJ::Sparticle> outgoing;
Event(
RHEJ::Sparticle in1, RHEJ::Sparticle in2,
std::initializer_list<RHEJ::Sparticle> out
):
incoming{std::move(in1), std::move(in2)},
outgoing{out}
{}
};
constexpr
RHEJ::pid::ParticleID PDG(int id){
return static_cast<RHEJ::pid::ParticleID>(id);
}
constexpr double alpha_s = 0.113559;
constexpr double mu = 125.;
constexpr double R = 0.4;
constexpr double min_jet_pt = 50.;
constexpr auto jet_def = fastjet::antikt_algorithm;
constexpr double ep = 1e-4;
void dump(Event const & ev){
LHEF::Writer writer{std::cerr};
std::cerr << std::setprecision(15);
RHEJ::UnclusteredEvent tmp;
tmp.incoming = ev.incoming;
tmp.outgoing = ev.outgoing;
tmp.central = {mu, mu, 0.};
RHEJ::Event out_ev{std::move(tmp), {jet_def, R}, min_jet_pt};
writer.hepeup = to_HEPEUP(std::move(out_ev), nullptr);
writer.writeEvent();
}
int main(){
std::vector<Event> events = {
+ /*
+ * reference matrix elements obtained from traditional HEJ (svn r3355)
+ * weights correspond to shat2*wt/(pdfPtr->pdfta*pdfPtr->pdftb)
+ * at the end of MEt2 in HJets.cxx
+ */
#include "ME_test_events.dat"
};
const std::vector<double> ref_weights{
1.20840366448943e-08,3.41482815561959e-18,4.09113550770465e-12,4.62704699562129e-15,4.48439126356004e-15,1.87972263060254e-12,1.44306589954128e-14,9.10660553817556e-13,6.78185181417658e-13,2.2313219279205e-14,5.04225136522914e-10,1.15445795400542e-13,4.25305992011228e-15,6.34464495221583e-11,8.1504882635896e-09,8.13050604628519e-10,1.16718768515638e-12,2.10641257506331e-08,6.62448031524987e-23,4.3319311109536e-07,1.20944494921349e-15,5.92539274973405e-05,4.4796169513938e-17,1.49273221119393e-15,2.43311569643869e-18,2.402863950606e-12,1.0620402696538e-11,5.35907540952987e-14,1.23031328706173e-12,2.79984692016006e-15,9.15681008462036e-19,6.82192590709917e-19,7.89237067022333e-13,1.86125787460036e-14,6.46979689455398e-07,3.04911830916725e-15,8.13163686285418e-15,1.53766124466997e-14,8.32276518653951e-14,8.88400099785712e-11,8.91102728936822e-06,1.80107644544759e-15,2.10582866365414e-09,3.10947796493188e-17,3.91243813845236e-15,3.26654480787734e-17,9.5447871679842e-14,2.59793669360488e-15,4.96134344054012e-13,5.81162371703576e-14,1.111167877034e-09,2.79109572058797e-16,4.46160513513727e-16,6.75218397576417e-15,2.68260561114305e-11,4.28989454505788e-16,5.8329868340234e-16,2.2024897389957e-09,2.57397955688743e-15,2.44654345522128e-14,2.30866752912176e-14,1.41827747056589e-07,1.08893147410797e-14,4.03382910318587e-11,3.1389869623674e-10,3.90505557430021e-13,1.08295237971538e-12,7.454591488958e-15,7.54483306433453e-14,1.02908259302562e-15,5.63531574394432e-14,6.30249429292698e-11,4.53762691808614e-13,1.89024041595359e-16,6.8158270448436e-15,1.14399844221906e-12,8.90256332512462e-11,6.06548092223796e-17,9.98362270779565e-12,4.59859716748231e-15,2.15099414771066e-10,2.11591109618363e-13,4.9742178019335e-09,1.33845681137089e-10,2.71186491627684e-16,8.86739836718398e-13,1.8249669781697e-14,1.92215499935624e-13,5.63538105112205e-14,4.14171568963499e-12,2.2717979957039e-15,5.55611156772101e-15,1.28364738716082e-16,5.3851134765245e-12,1.24211103260246e-08,4.83232530709888e-17,0.000521569751168655,1.6806219324588e-15,6.59977348355555e-11,9.1627970683281e-14,2.76820118228537e-12,4.0198982918336e-11,1.48426883756103e-18,1.40114232815435e-18,1.7965991343524e-11,7.6977628175923e-10,4.62196685302415e-14,2.69246069544232e-08,1.51516950568948e-13,6.99496684443045e-13,5.18573546306991e-16,0.0434630460883225,1.9837330600509e-23,2.67017601462925e-12,9.8878856049861e-16,8.72580364190897e-11,9.78094461016429e-06,8.71781003862465e-17,1.1559862791241e-08,4.84083438123941e-11,0.0139742563571754,1.99847601444777e-13,3.04915596353023e-17,5.93151845056274e-16,2.17435618150856e-14,5.33225221812805e-08,3.87773807827851e-14,4.50810202343348e-08,1.21607462589437e-08,3.1594598104258e-13,3.28656273761685e-16,2.50210924939855e-21,2.38126296013289e-13,1.25433376653286e-11,2.14322725557091e-14,1.49146107213853e-09,1.242585653879e-15,3.55954759815053e-17,3.58323190858796e-10,4.51037144634267e-14,7.55931651143187e-09,5.37698824199559e-14,7.38868781004529e-14,1.69357650230295e-11,5.80639629823124e-10,7.4837350747944e-17,6.96956083393151e-16,1.81834207579059e-10,2.82540028744915e-12,4.0395726561383e-13,1.79694539029015e-12,1.31041890677814e-16,1.30438189219159e-07,1.43221371664219e-10,9.114047729275e-16,7.34577355123687e-09,3.22593068239571e-15,4.15681922353836e-15,3.36408379142654e-16,1.44480854726148e-11,1.01363611458694e-13,1.00219358430676e-12,0.00338137229726986,4.29548017717404e-15,2.60395354258443e-14,4.1476296627427e-14,1.06855911096978e-13,8.22439039784712e-11,2.73629652512136e-10,5.15659349602634e-11,1.01584718933351e-13,6.94380398977803e-16,2.18931000347143e-18,3.07933496617785e-09,1.66390942837878e-16,3.52525907142114e-13,3.78080647432571e-14,1.17556440172277e-12,6.59751630435329e-17,1.21755902521401e-08,2.55450721249666e-11,1.53256550428156e-14,1.08380893991568e-12,4.16098986326189e-13,2.90205083341855e-14,5.511460594056e-18,1.57699398911427e-16,7.00104401760928e-13,3.58046798793545e-12,2.01229361268661e-12,3.73614518690553e-07,2.14710426845508e-14,3.70105600162546e-13,7.98071394035332e-15,7.37834963791502e-16,1.60529982053429e-09,1.08256973923745e-12,1.0334353499219e-13,9.49543099778791e-11,1.9115591364306e-12,4.0448657562838e-16,8.71422438179998e-15,1.25631778664749e-16,3.14479096095934e-19,2.64428535684163e-16,1.12107615349188e-11,1.23625298812219e-13,4.01428912903254e-13,2.70797308155832e-11,6.56279497820654e-09,7.13734059806893e-16,0.000782868291989126,3.76614026148548e-16,1.26257000484765e-13,9.58631950869197e-13,8.31708875806124e-11,2.34133301418612e-13,9.56863619261736e-14,4.64712032868663e-12,2.86185199150703e-12,1.33576465609689e-11,1.10079981208014e-15,1.0127254462106e-12,1.1462054990501e-12,6.6945137657857e-11,5.58956247731427e-17,9.55921087729154e-12,2.34217482968982e-13,9.06915368555449e-19,5.71005344841003e-15,7.30755001164554e-09,6.21879071656432e-07,3.18226338142145e-17,7.85187718829052e-15,1.12784019414768e-13,1.26620458266236e-12,3.69177699852588e-14,5.15075992359208e-08,4.92072565553936e-15,1.85781482577606e-18,4.45444190278706e-16,1.26266755594178e-13,2.4288203970545e-16,3.05176309462854e-13,1.10490541122231e-16,8.47749453725764e-16,5.15877353528789e-11,5.24254149207253e-17,1.91355494316452e-10,5.00556344079427e-17,2.17154626918891e-15,8.84358943960976e-15,3.61414423657153e-11,4.82735747825927e-09,1.29064865560832e-13,8.42127739585314e-10,3.92906216459766e-16,6.04446459041361e-17,5.89524177334148e-12,5.46194519536606e-19,1.34088171433487e-11,5.52860041827642e-15,1.75477788501097e-15,1.73667838135065e-16,1.65397026290053e-13,1.12185521855334e-09,6.08386662697057e-17,4.59141650497559e-07,1.5209325219688e-11,2.23037744763897e-13,2.01952877105987e-13,6.8597135743949e-17,1.21998302502671e-11,3.21171910404898e-17,5.26021980482437e-15,7.19705632175135e-19,1.20963419308456e-09,6.88528237531901e-13,1.98848279147693e-11,4.01262000581642e-12,5.20247739003617e-13,1.30428205622828e-14,4.11298291694491e-13,5.52043003142038e-10,1.23093415613032e-16,7.36616975809533e-17,3.45364902177034e-16,4.34985181424009e-08,5.38804276057675e-14,9.48300760869968e-13,1.22337104642822e-09,8.35787040119371e-14,1.55423397182344e-15,1.42203125672685e-13,3.18903277676011e-13,3.88652741986258e-15,7.04896809011702e-18,1.32646740227354e-21,7.99231856556843e-13,5.95802234870006e-09,1.65149306543534e-12,1.33064116256099e-11,1.56832702886549e-12,1.10299934631458e-08,3.67898989067657e-16,7.84682358054423e-17,2.57971628718788e-09,8.33182695406546e-17,1.18568992167241e-10,1.84355449893326e-12,2.55066898328171e-13,2.83673294248314e-14,1.78606190612042e-12,3.58348766614931e-12,4.59943073492537e-16,9.51522624162865e-13,8.0286195911246e-14,1.36695977212813e-16,2.66996900743536e-13,1.18874419459462e-15,4.20697505402443e-07,1.44157427471053e-12,3.88352207261609e-12,1.10352682706254e-15,9.16382540998255e-16,1.6918487768822e-14,6.27001020277801e-13,1.19982339734628e-15,2.08524585760556e-15,3.98635216491517e-11,2.29164410091261e-16,6.59993812570474e-14,9.09720299660866e-19,1.00791143935241e-13,3.48282488442434e-15,3.82462701684145e-11,1.82912920125976e-14,2.35934504263123e-12,3.07876467756239e-11,3.01958860729077e-05,1.34442938440306e-12,1.98172170171169e-11,2.54803366934403e-14,4.63525528427102e-13,1.98623090357032e-11,1.17532063246858e-12,1.62645374532443e-16,2.23093836031326e-11,5.23718977936426e-13,6.30254161980384e-15,5.64457070654369e-16,5.8598582297477e-15,1.54121478101166e-11,2.54757492800786e-11,7.24784320481052e-13,7.10590291356503e-11,1.90012029508776e-12,5.38421849395334e-05,9.2027530054952e-13,0.00802523961864568,1.69971085315301e-12,4.17128853468654e-15,2.58531014487319e-15,6.72871249369096e-12,3.68392430578061e-16,1.56192266889718e-14,8.09758960836169e-16,1.33670176619002e-07,3.40402526854793e-10,7.14029639915786e-15,1.15907463202726e-14
};
const std::vector<double> ref_virt{
0.000407437434887122,0.0616303680115081,0.0619338272462471,0.00113600815259533,0.0611805205486735,0.000565698082480416,0.516976611025519,0.0361140019339901,0.00145895360905874,0.134071980986426,0.000248975936203567,0.199294625870235,0.24358305091454,4.74361172738652e-05,5.89666196640687e-05,0.00531027816302915,0.000200243383490462,0.0805261711344379,0.00120027603203917,0.00201705178179895,0.339181495419129,0.000229863490455131,0.0523929260778993,0.000263057002724373,0.00244241108460007,0.0261467924657004,5.9657372131903e-05,0.200512468910312,0.00810099141722299,0.013973033333956,0.0108154079325531,0.158574821336736,0.00070998197888614,0.000498608284914822,3.73054997210496e-05,0.0589374030166269,0.00230622721548646,0.0339477696166811,0.000912463632412847,0.000736710687806361,0.00409015383244036,0.245951477936874,0.0837392661445892,0.222078993106222,0.00532590926455614,0.625630979237125,0.0496003388518398,0.10116468200925,0.0661786669164834,0.0210544045965245,0.0244625178289293,0.0146510511459362,0.101257018174142,0.00316763508016474,0.00143243489499069,0.333592621341247,0.00370907687545401,0.00213872056006003,0.01900030889465,0.167144065343861,2.6100961817368e-05,4.91550277446228e-05,0.0166499896061438,0.00011398434602661,0.000359768751250008,0.00225147558163185,0.0317147645964709,0.0842595514695342,0.35407612826878,0.00010526694382132,0.0430608610020923,3.5715013356257e-05,0.000838215519554185,0.00192018265054966,0.0713712284182796,0.000485232329561933,0.00152471141515176,0.000311862603200811,0.166619420201183,0.0199021050517686,1.41779057397457e-05,0.0017318628596011,0.00053760812759845,0.0686125677803084,0.0532620439395034,0.00166666380089786,0.0113533421860016,0.000109349820296255,0.0115667016173654,0.00173679164546385,0.0256607191277514,0.00120155627964714,0.0111815333349668,0.00349052851306921,0.00531735228893806,0.00152975285273546,4.66037248418355e-05,0.0425850611539657,0.00117697776465851,1.23387751381228e-05,0.0038237470220905,0.0905644410585002,0.02820127090728,0.00024367136851271,0.147636957219027,0.000593136877076053,0.0618471045067722,9.07687985864349e-06,0.000158960279301719,0.0887358461047651,0.0332353914977947,7.4341716367417e-06,0.000892008177580303,0.000741022113249675,0.00237741099974918,0.00354813689304541,0.00074173961384706,0.00467141577324516,0.00842296482521057,0.00821661392007147,4.83792120557372e-05,0.00425252516494632,0.000432204419963046,0.00586736518573752,0.000956530707046858,0.000236851602228008,0.0892427215701044,1.38830962250345e-05,0.0080267826500828,0.284615916615109,0.000228095910528196,0.000135342214864163,7.55341699407512e-06,0.000571884401074523,0.0101768016413147,0.0171951107578717,0.00189669876455143,0.00394930354845142,0.00328294862282636,0.00123193340562429,0.00331783483858544,0.00319051646342221,0.00155129659005506,8.87358953340508e-05,0.27623170363077,0.000223827470760779,0.0568015505265356,0.000286151801691539,0.0193226057191209,0.00123755083758374,0.000465248541258071,0.0760142560624272,0.00112490607994365,0.000144842352216183,0.197405681090406,0.00129883979140284,0.0930965339402953,0.0010308437485118,0.670013359364299,0.00432667098624779,0.0454915871138651,0.00251703872973342,0.00423441033466031,0.000345157767746845,0.00161415356368555,0.00416069958683314,0.00164939718639514,0.0111406419016027,0.00655421252542946,0.00851532238518754,0.00348709627275328,0.496059905302338,0.0371917367183897,0.000211531938236732,0.000258942974591912,0.17759337052515,0.00418622722335427,0.166120697713642,0.00360502919950296,0.000407414278464457,0.0360738279541873,0.0288893174873495,0.00047259015965458,1.56547248986137e-05,0.834956154090444,0.00732703650642937,0.00474657800944094,0.00273096554923287,0.257207154892208,0.00469808491585188,0.000270576144489275,0.000136075839471612,0.097205919234904,0.000607126031374257,0.000136619984636834,0.000160087925098193,0.0073858227644151,0.538982397271589,0.00336284303209108,0.00152656630103797,0.367113302445716,7.48038626150403e-05,2.60037367946894e-05,4.68016275421979e-05,0.000158257654442151,0.000139074963089943,0.000665698879408135,0.000343142682314476,0.0686327476602143,0.000242720484154181,0.0343708294206813,0.00114889968298977,0.0586454267607763,0.0681515891243611,0.00825429790468652,0.023992046337532,0.0176721744671367,0.000118812692493978,0.152442537631382,0.017153140472422,0.0064172004743378,3.15311961394686e-05,0.0176719977165195,0.00523267370956132,0.00696013339412051,0.0283129688131134,0.000387153677372906,0.00691199037138631,0.000146275451696295,0.701230156661611,0.000192846626056869,0.000164329434653947,0.00161262788918379,0.138649850457724,0.0298984433917698,0.0318122273046375,3.04827880239274e-05,0.000447390422746664,0.00113026494701542,0.00113285595566244,2.89098339981399e-05,0.00113856357833898,0.0539144571449436,0.0708152617112735,0.000146823107858066,0.0226524231026193,0.00337461891944858,0.000951319767163979,0.00173528759879105,0.027553837556441,0.000644790987241654,0.00202675597720372,0.0904812026317863,0.000121859990107361,0.0778263959001868,0.0078676495423855,0.00106041104696191,0.0460181415526085,0.111418229601818,0.00485578450515022,0.000301615731939099,0.0991074788208953,0.000545043001929007,0.00699096247420379,0.203553031240112,0.00129391506623192,0.0163036481873179,9.84784360850182e-06,0.0011930992876943,0.000233138513444116,0.00327134762224621,0.00934510369593952,0.0257690022271133,0.0586620958023085,0.0044269686967688,0.00458204111385776,0.000575173870997244,0.0196091953126155,0.0172446254918132,0.0023874923763825,0.00203444177442346,0.00590189106361814,0.160798676736245,0.000218465854877618,0.0449130791723244,0.0221907455403962,0.00577930006864582,4.68129848124974e-05,4.39088532565896e-05,0.0957291068496414,0.000222196152405895,0.00244492299475677,0.0649283175892222,0.00246575641284538,0.00395895850035338,0.452609135498808,2.18247701431579e-05,0.000103691495491759,0.0300254483476907,0.000997825085990438,0.00120168803041095,0.000488310402876861,5.02375648302248e-05,0.0889218155230037,0.886196353445691,0.537472442381394,0.000989158503397272,0.00585452117968649,0.000166646041974894,0.0559913557856292,0.271624614149917,1.9604193016588e-05,0.00400455993655262,0.000719953658022318,0.00769221521719457,0.0122222590022294,0.00869539856191144,0.00665334452394116,0.00269456106141176,0.0203144580365219,1.06269256412071e-05,0.00266305904353431,0.000866750257077464,0.0118678506617259,0.00774140571086026,0.472668086097713,0.00449314655439557,0.000363054262255041,0.000123958362353204,0.0593128481039698,6.49758497533422e-05,0.0147733538718899,0.000181422709422017,0.124358363237549,0.0149418567595413,0.000496790567694211,0.00319609174497007,0.171371757877329,0.000132314177013318,0.000118421517677613,0.0384411149022143,0.0560023977555578,0.000545951922843868,0.00923849055867683,7.42750022847054e-05,0.0034657141222543,0.0814951536826408,0.000247436489650107,0.00279135741337402,0.128217103654628,0.000239656955359697,0.0286081920788208,0.0264861644616376,0.000700223039958121,0.000224261722502375,0.206999869456071,0.0245477336638465,0.00020080171595025,0.000534370708560753,0.00114853261836608,0.0140588914482069,0.0491695416950957,0.0219392429988939,0.000187573794241001,2.76073824918907e-05,0.0177905351398601,0.0227102533101537,0.0156338659659345,0.0045130866637909,0.0118789100672146,0.0684121072806861
};
RHEJ::MatrixElement ME{
+ [](double){ return alpha_s; },
{jet_def, R}, min_jet_pt, false,
RHEJ::HiggsCouplingSettings{}
};
assert(events.size() == ref_weights.size());
for(size_t i = 0; i < events.size(); ++i){
std::sort(
begin(events[i].outgoing), end(events[i].outgoing),
RHEJ::rapidity_less{}
);
const double our_ME = ME.tree(
- alpha_s, mu, events[i].incoming, events[i].outgoing, true
+ mu, events[i].incoming, events[i].outgoing, true
);
const double deviation = our_ME/ref_weights[i] - 1.;
if(std::abs(deviation) > ep){
std::cerr
<< "Matrix element deviates by factor of " << 1+deviation << ":\n"
"Is " << our_ME << ", should be " << ref_weights[i] << "\n"
"Event:\n";
dump(events[i]);
return EXIT_FAILURE;
}
const double our_virt = ME.virtual_corrections(
- alpha_s, mu, events[i].incoming, events[i].outgoing
+ mu, events[i].incoming, events[i].outgoing
);
const double virt_deviation = our_virt/ref_virt[i] - 1.;
if(std::abs(virt_deviation) > ep){
std::cerr
<< "Virtual corrections deviate by factor of " << 1+virt_deviation << ":\n"
"Is " << our_virt << ", should be " << ref_virt[i] << "\n"
"Event:\n";
dump(events[i]);
return EXIT_FAILURE;
}
}
}
diff --git a/t/test_ME_hjets_mt174.cc b/t/test_ME_hjets_mt174.cc
index 0aaa81b..e46afbc 100644
--- a/t/test_ME_hjets_mt174.cc
+++ b/t/test_ME_hjets_mt174.cc
@@ -1,78 +1,84 @@
#include "LHEF/LHEF.h"
#include "RHEJ/MatrixElement.hh"
#include "RHEJ/Event.hh"
struct Event{
double weight;
std::array<RHEJ::Sparticle, 2> incoming;
std::vector<RHEJ::Sparticle> outgoing;
Event(
double wt,
RHEJ::Sparticle in1, RHEJ::Sparticle in2,
std::initializer_list<RHEJ::Sparticle> out
):
weight{wt},
incoming{std::move(in1), std::move(in2)},
outgoing{out}
{}
};
constexpr
RHEJ::pid::ParticleID PDG(int id){
return static_cast<RHEJ::pid::ParticleID>(id);
}
constexpr double alpha_s = 0.112654;
constexpr double mu = 125.;
constexpr double R = 0.4;
constexpr double min_jet_pt = 35.;
constexpr auto jet_def = fastjet::antikt_algorithm;
constexpr double ep = 1e-4;
constexpr double alpha_s_mH_RHEJ = 0.113559;
void dump(Event const & ev){
LHEF::Writer writer{std::cerr};
std::cerr << std::setprecision(15);
RHEJ::UnclusteredEvent tmp;
tmp.incoming = ev.incoming;
tmp.outgoing = ev.outgoing;
tmp.central = {mu, mu, 0.};
RHEJ::Event out_ev{std::move(tmp), {jet_def, R}, min_jet_pt};
writer.hepeup = to_HEPEUP(std::move(out_ev), nullptr);
writer.writeEvent();
}
int main(){
std::vector<Event> events = {
+ /*
+ * reference matrix elements obtained from traditional HEJ (svn r3355)
+ * weights correspond to shat2*wt/(pdfPtr->pdfta*pdfPtr->pdftb)
+ * at the end of MEt2 in HJets.cxx
+ */
#include "ME_hjets_mt174.dat"
};
RHEJ::HiggsCouplingSettings settings;
settings.mt = 174.;
settings.use_impact_factors = false;
settings.mb = 5.;
settings.include_bottom = true;
RHEJ::MatrixElement ME{
+ [](double){ return alpha_s; },
{jet_def, R}, min_jet_pt, false, settings
};
for(size_t i = 0; i < events.size(); ++i){
std::sort(
begin(events[i].outgoing), end(events[i].outgoing),
RHEJ::rapidity_less{}
);
const double our_ME = ME.tree(
- alpha_s, mu, events[i].incoming, events[i].outgoing, true
- )*alpha_s*alpha_s/(alpha_s_mH_RHEJ*alpha_s_mH_RHEJ);
+ mu, events[i].incoming, events[i].outgoing, true
+ );
const double deviation = our_ME/events[i].weight - 1.;
if(std::abs(deviation) > ep){
std::cerr
<< "Matrix element deviates by factor of " << 1+deviation << ":\n"
"Is " << our_ME << ", should be " << events[i].weight << "\n"
"Event:\n";
dump(events[i]);
return EXIT_FAILURE;
}
}
}

File Metadata

Mime Type
text/x-diff
Expires
Tue, Jan 21, 2:18 AM (1 d, 21 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
4243644
Default Alt Text
(288 KB)

Event Timeline