diff --git a/CMakeLists.txt b/CMakeLists.txt index 684d4da..8de4ef5 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -1,376 +1,384 @@ cmake_minimum_required(VERSION 3.1 FATAL_ERROR) set(CMAKE_LEGACY_CYGWIN_WIN32 0) set(CMAKE_EXPORT_COMPILE_COMMANDS ON) project("HEJ" VERSION 2.0.5 LANGUAGES C CXX) # Set a default build type if none was specified set(default_build_type "RelWithDebInfo") if(NOT CMAKE_BUILD_TYPE AND NOT CMAKE_CONFIGURATION_TYPES) message(STATUS "Setting build type to '${default_build_type}' as none was specified.") set(CMAKE_BUILD_TYPE "${default_build_type}" CACHE STRING "Choose the type of build." FORCE) # Set the possible values of build type for cmake-gui set_property(CACHE CMAKE_BUILD_TYPE PROPERTY STRINGS "Debug" "Release" "MinSizeRel" "RelWithDebInfo") endif() ## Flags for the compiler. No warning allowed. if (CMAKE_COMPILER_IS_GNUCC OR CMAKE_COMPILER_IS_GNUCXX) set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wall -Wextra") elseif (MSVC) set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} /W4 /WX /EHsc") endif() set(CMAKE_CXX_STANDARD_REQUIRED ON) set(CMAKE_CXX_STANDARD 14) ## Create Version set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} "${CMAKE_SOURCE_DIR}/cmake/Modules/") # Get the latest abbreviated commit hash of the working branch execute_process( COMMAND git rev-parse HEAD WORKING_DIRECTORY ${CMAKE_CURRENT_SOURCE_DIR} OUTPUT_VARIABLE PROJECT_GIT_REVISION OUTPUT_STRIP_TRAILING_WHITESPACE ) # Get the current working branch execute_process( COMMAND git rev-parse --abbrev-ref HEAD WORKING_DIRECTORY ${CMAKE_CURRENT_SOURCE_DIR} OUTPUT_VARIABLE PROJECT_GIT_BRANCH OUTPUT_STRIP_TRAILING_WHITESPACE ) ## target directories for install set(INSTALL_INCLUDE_DIR_BASE include) set(INSTALL_INCLUDE_DIR ${INSTALL_INCLUDE_DIR_BASE}/HEJ) set(INSTALL_BIN_DIR bin) set(INSTALL_LIB_DIR lib) set(INSTALL_CONFIG_DIR lib/cmake/HEJ) ## Template files configure_file( ${CMAKE_CURRENT_SOURCE_DIR}/cmake/Templates/Version.hh.in ${PROJECT_BINARY_DIR}/include/HEJ/Version.hh @ONLY ) configure_file( ${CMAKE_CURRENT_SOURCE_DIR}/cmake/Templates/HEJ-config.cc.in ${PROJECT_BINARY_DIR}/src/bin/HEJ-config.cc @ONLY ) # Generate CMake config file include(CMakePackageConfigHelpers) configure_package_config_file( cmake/Templates/hej-config.cmake.in ${CMAKE_CURRENT_BINARY_DIR}/${INSTALL_CONFIG_DIR}/hej-config.cmake INSTALL_DESTINATION ${INSTALL_CONFIG_DIR} PATH_VARS INSTALL_INCLUDE_DIR_BASE INSTALL_LIB_DIR ) write_basic_package_version_file( ${CMAKE_CURRENT_BINARY_DIR}/${INSTALL_CONFIG_DIR}/hej-config-version.cmake COMPATIBILITY SameMajorVersion ) install( FILES ${CMAKE_CURRENT_BINARY_DIR}/${INSTALL_CONFIG_DIR}/hej-config.cmake ${CMAKE_CURRENT_BINARY_DIR}/${INSTALL_CONFIG_DIR}/hej-config-version.cmake DESTINATION ${INSTALL_CONFIG_DIR}) ## Add directories and find dependences include_directories(${CMAKE_CURRENT_SOURCE_DIR}/include ${PROJECT_BINARY_DIR}/include) find_package(fastjet REQUIRED) include_directories(${fastjet_INCLUDE_DIRS}) find_package(CLHEP 2.3 REQUIRED) include_directories(${CLHEP_INCLUDE_DIRS}) find_package(LHAPDF REQUIRED) include_directories(${LHAPDF_INCLUDE_DIRS}) ## Amend unintuitive behaviour of FindBoost.cmake ## Priority of BOOST_ROOT over BOOSTROOT matches FindBoost.cmake ## at least for cmake 3.12 if(DEFINED BOOST_ROOT) message("BOOST_ROOT set - only looking for Boost in ${BOOST_ROOT}") set(Boost_NO_BOOST_CMAKE ON) elseif(DEFINED BOOSTROOT) message("BOOSTROOT set - only looking for Boost in ${BOOSTROOT}") set(Boost_NO_BOOST_CMAKE ON) endif() find_package(Boost REQUIRED COMPONENTS iostreams) include_directories(${Boost_INCLUDE_DIRS}) find_package(yaml-cpp) # requiring yaml does not work with fedora include_directories(${YAML_CPP_INCLUDE_DIR}) if(${EXCLUDE_HepMC}) message(STATUS "Skipping HepMC") # avoid "unused variable" warning if EXCLUDE_rivet is set by user set(EXCLUDE_rivet TRUE) else() find_package(HepMC 2) endif() if(${HepMC_FOUND}) set( CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -DHEJ_BUILD_WITH_HepMC_VERSION=${HepMC_VERSION_MAJOR}" ) include_directories(${HepMC_INCLUDE_DIRS}) if(${EXCLUDE_rivet}) message(STATUS "Skipping rivet") else() find_package(rivet) endif() if(${rivet_FOUND}) include_directories(${rivet_INCLUDE_DIRS}) set( CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -DHEJ_BUILD_WITH_RIVET" ) endif() endif() if(${EXCLUDE_QCDloop}) message(STATUS "Skipping QCDloop") else() find_package(QCDloop 2) endif() if(${QCDloop_FOUND}) set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -DHEJ_BUILD_WITH_QCDLOOP") include_directories(SYSTEM ${QCDloop_INCLUDE_DIRS}) endif() if(${EXCLUDE_HighFive}) message(STATUS "Skipping HighFive") else() find_package(HighFive QUIET) find_package_handle_standard_args( HighFive CONFIG_MODE ) endif() if(${HighFive_FOUND}) set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -DHEJ_BUILD_WITH_HDF5") include_directories($<TARGET_PROPERTY:HighFive,INTERFACE_INCLUDE_DIRECTORIES>) endif() add_subdirectory(src) ## define executable add_executable(HEJ src/bin/HEJ.cc) ## link libraries target_link_libraries(HEJ hejlib) add_executable(HEJ-config src/bin/HEJ-config.cc) file(GLOB hej_headers ${CMAKE_CURRENT_SOURCE_DIR}/include/HEJ/*.hh ${PROJECT_BINARY_DIR}/include/HEJ/*.hh ) file(GLOB hej_detail_headers ${CMAKE_CURRENT_SOURCE_DIR}/include/HEJ/detail/*.hh) file(GLOB lhef_headers ${CMAKE_CURRENT_SOURCE_DIR}/include/LHEF/*.h) install(FILES ${hej_headers} DESTINATION ${INSTALL_INCLUDE_DIR}) install(FILES ${hej_detail_headers} DESTINATION ${INSTALL_INCLUDE_DIR}/detail/) install(FILES ${lhef_headers} DESTINATION ${INSTALL_INCLUDE_DIR_BASE}/LHEF/) install(TARGETS HEJ HEJ-config DESTINATION ${INSTALL_BIN_DIR}) ## tests enable_testing() set(tst_dir "${CMAKE_CURRENT_SOURCE_DIR}/t") set(tst_ME_data_dir "${tst_dir}/ME_data") # test event classification add_executable(test_classify ${tst_dir}/test_classify.cc) target_link_libraries(test_classify hejlib) add_test( NAME t_classify COMMAND test_classify ${tst_dir}/classify.lhe.gz ) add_executable(test_classify_new ${tst_dir}/test_classify_new.cc) target_link_libraries(test_classify_new hejlib) add_test( NAME t_classify_new COMMAND test_classify_new ) # test phase space point add_executable(test_psp ${tst_dir}/test_psp.cc) target_link_libraries(test_psp hejlib) add_test( NAME t_psp COMMAND test_psp ${tst_dir}/psp_gen.lhe.gz ) # test importing scales add_library(scales SHARED ${tst_dir}/scales.cc) add_executable(test_scale_import ${tst_dir}/test_scale_import) target_link_libraries(test_scale_import hejlib) add_test( NAME t_scale_import COMMAND test_scale_import ${tst_dir}/jet_config_with_import.yml ) # test scale arithmetic (e.g. 2*H_T/4) add_executable(test_scale_arithmetics ${tst_dir}/test_scale_arithmetics) target_link_libraries(test_scale_arithmetics hejlib) add_test( NAME t_scale_arithmetics COMMAND test_scale_arithmetics ${tst_dir}/jet_config.yml ${tst_dir}/2j.lhe.gz ) # test "ParameterDescription" add_executable(test_descriptions ${tst_dir}/test_descriptions) target_link_libraries(test_descriptions hejlib) add_test( NAME t_descriptions COMMAND test_descriptions ) # test "EventParameters*Weight" add_executable(test_parameters ${tst_dir}/test_parameters) target_link_libraries(test_parameters hejlib) add_test( NAME test_parameters COMMAND test_parameters ) # test colour generation add_executable(test_colours ${tst_dir}/test_colours) target_link_libraries(test_colours hejlib) add_test( NAME t_colour_flow COMMAND test_colours ) # test matrix elements add_executable(test_ME_generic ${tst_dir}/test_ME_generic.cc) target_link_libraries(test_ME_generic hejlib) add_test( NAME t_ME_j COMMAND test_ME_generic ${tst_ME_data_dir}/config_mtinf.yml ${tst_ME_data_dir}/ME_jets.dat ${tst_ME_data_dir}/PSP_jets.lhe.gz ) add_test( NAME t_ME_h COMMAND test_ME_generic ${tst_ME_data_dir}/config_mtinf.yml ${tst_ME_data_dir}/ME_h_mtinf.dat ${tst_ME_data_dir}/PSP_h.lhe.gz ) if(${QCDloop_FOUND}) add_test( NAME t_ME_h_mt COMMAND test_ME_generic ${tst_ME_data_dir}/config_mt.yml ${tst_ME_data_dir}/ME_h_mt_tree.dat ${tst_ME_data_dir}/PSP_h.lhe.gz ) add_test( NAME t_ME_h_mtmb COMMAND test_ME_generic ${tst_ME_data_dir}/config_mtmb.yml ${tst_ME_data_dir}/ME_h_mtmb_tree.dat ${tst_ME_data_dir}/PSP_h.lhe.gz ) endif() add_test( NAME t_ME_w_FKL COMMAND test_ME_generic ${tst_ME_data_dir}/config_w_ME.yml ${tst_ME_data_dir}/ME_w_FKL_tree.dat ${tst_ME_data_dir}/PSP_w_FKL.lhe.gz ) add_test( NAME t_ME_w_FKL_virt COMMAND test_ME_generic ${tst_ME_data_dir}/config_w_ME.yml ${tst_ME_data_dir}/ME_w_FKL_virt.dat ${tst_ME_data_dir}/PSP_w_FKL.lhe.gz ) add_test( NAME t_ME_Wp COMMAND test_ME_generic ${tst_ME_data_dir}/config_w_ME.yml ${tst_ME_data_dir}/ME_Wp.dat ${tst_ME_data_dir}/PSP_Wp.lhe.gz ) add_test( NAME t_ME_Wm COMMAND test_ME_generic ${tst_ME_data_dir}/config_w_ME.yml ${tst_ME_data_dir}/ME_Wm.dat ${tst_ME_data_dir}/PSP_Wm.lhe.gz ) # test main executable if(${HepMC_FOUND}) file(READ "${tst_dir}/jet_config.yml" config) file(WRITE "${tst_dir}/jet_config_withHepMC.yml" "${config} - tst.hepmc") if(${rivet_FOUND}) file(READ "${tst_dir}/jet_config_withHepMC.yml" config) file(WRITE "${tst_dir}/jet_config_withRivet.yml" "${config}\n\nanalysis:\n rivet: MC_XS\n output: tst") add_test( NAME t_main COMMAND HEJ ${tst_dir}/jet_config_withRivet.yml ${tst_dir}/2j.lhe.gz ) else() # no rivet add_test( NAME t_main COMMAND HEJ ${tst_dir}/jet_config_withHepMC.yml ${tst_dir}/2j.lhe.gz ) endif() if(${HepMC_VERSION_MAJOR} GREATER 2) add_executable(check_hepmc ${tst_dir}/check_hepmc.cc) target_link_libraries(check_hepmc hejlib) # check that HepMC output is correct add_test( NAME t_hepmc COMMAND check_hepmc tst.hepmc ) endif() else() # no HepMC add_test( NAME t_main COMMAND HEJ ${tst_dir}/jet_config.yml ${tst_dir}/2j.lhe.gz ) endif() # check that LHEF output is correct add_executable(check_lhe ${tst_dir}/check_lhe.cc) target_link_libraries(check_lhe hejlib) add_test( NAME t_lhe COMMAND check_lhe tst.lhe ) # check HDF5 reader if(${HighFive_FOUND}) add_executable(test_hdf5 ${tst_dir}/test_hdf5.cc) target_link_libraries(test_hdf5 hejlib) add_test( NAME t_hdf5 COMMAND test_hdf5 ${tst_dir}/Wm9-g4-repack.hdf5 ) endif() # test boson reconstruction add_executable(cmp_events ${tst_dir}/cmp_events.cc) target_link_libraries(cmp_events hejlib) add_test( NAME t_epnu_2j_noW COMMAND cmp_events ${tst_dir}/epnu2jLOFKL_unweight.lhe.tar.gz ${tst_dir}/epnu2jLOFKL_unweight_noW.lhe.gz ) # test resummed result add_executable(check_res ${tst_dir}/check_res.cc) target_link_libraries(check_res hejlib) if(${TEST_ALL}) # deactivate long tests by default add_test( NAME t_2j COMMAND check_res ${tst_dir}/2j.lhe.gz 3.49391e+07 419684 ) add_test( NAME t_3j COMMAND check_res ${tst_dir}/3j.lhe.gz 2.37902e+06 25746.6 ) add_test( NAME t_4j COMMAND check_res ${tst_dir}/4j.lhe.gz 603713 72822.6 ) add_test( NAME t_h_3j COMMAND check_res ${tst_dir}/h_3j.lhe.gz 0.821622 0.0220334 ) add_test( NAME t_h_3j_uno COMMAND check_res ${tst_dir}/h_3j_uno.lhe.gz 0.0261968 0.000341549 uno ) add_test( NAME t_epnu_2j COMMAND check_res ${tst_dir}/epnu2jLOFKL_unweight.lhe.tar.gz 262.7 3 ) add_test( NAME t_MGepnu_3j COMMAND check_res ${tst_dir}/MGepnu3j_unweighted.tar.gz 38.9512 1 ) add_test( NAME t_MGemnubar_3j COMMAND check_res ${tst_dir}/MGemnubar3j_unweighted.tar.gz 24.1575 1 ) add_test( + NAME t_MGepnu_3j_split + COMMAND check_res ${tst_dir}/MGepnu3j_unweighted.tar.gz 6.51853 0.1 split + ) + add_test( NAME t_MGepnu_4j COMMAND check_res ${tst_dir}/MGepnu4j_unweighted.tar.gz 10.2542 0.135106 ) add_test( NAME t_MGemnubar_4j COMMAND check_res ${tst_dir}/MGemnubar4j_unweighted.tar.gz 5.57909 0.0300496 ) + add_test( + NAME t_MGemnubar_4j_qqxmid + COMMAND check_res ${tst_dir}/MGemnubar4j_unweighted.tar.gz 0.732084 0.005 qqxmid + ) endif() diff --git a/t/check_res.cc b/t/check_res.cc index 8e26730..4a97302 100644 --- a/t/check_res.cc +++ b/t/check_res.cc @@ -1,141 +1,152 @@ /** * \authors The HEJ collaboration (see AUTHORS for details) * \date 2019 * \copyright GPLv2 or later */ #include <iostream> #include "LHEF/LHEF.h" #include "HEJ/Event.hh" #include "HEJ/EventReweighter.hh" #include "HEJ/Mixmax.hh" #include "HEJ/stream.hh" #define ASSERT(x) if(!(x)) { \ std::cerr << "Assertion '" #x "' failed.\n"; \ return EXIT_FAILURE; \ } namespace{ const fastjet::JetDefinition jet_def{fastjet::kt_algorithm, 0.4}; const fastjet::JetDefinition Born_jet_def{jet_def}; constexpr double Born_jetptmin = 30; constexpr double extpartonptmin = 30; constexpr double max_ext_soft_pt_fraction = std::numeric_limits<double>::infinity(); constexpr double jetptmin = 35; constexpr bool log_corr = false; using EventTreatment = HEJ::EventTreatment; using namespace HEJ::event_type; HEJ::EventTreatMap treat{ {no_2_jets, EventTreatment::discard}, {bad_final_state, EventTreatment::discard}, {FixedOrder, EventTreatment::discard}, {unof, EventTreatment::discard}, {unob, EventTreatment::discard}, {qqxexb, EventTreatment::discard}, {qqxexf, EventTreatment::discard}, {qqxmid, EventTreatment::discard}, {FKL, EventTreatment::reweight} }; /// true if colour is allowed for particle bool correct_colour(HEJ::Particle const & part){ if(HEJ::is_AWZH_boson(part) && !part.colour) return true; if(!part.colour) return false; int const colour = part.colour->first; int const anti_colour = part.colour->second; if(part.type == HEJ::ParticleID::gluon) return colour != anti_colour && colour > 0 && anti_colour > 0; if(HEJ::is_quark(part)) return anti_colour == 0 && colour > 0; return colour == 0 && anti_colour > 0; } bool correct_colour(HEJ::Event const & ev){ if(!HEJ::event_type::is_HEJ(ev.type())) return true; for(auto const & part: ev.incoming()){ if(!correct_colour(part)) return false; } for(auto const & part: ev.outgoing()){ if(!correct_colour(part)) return false; } return true; } }; int main(int argn, char** argv) { if(argn == 5 && std::string(argv[4]) == "uno"){ --argn; treat[unof] = EventTreatment::reweight; treat[unob] = EventTreatment::reweight; treat[FKL] = EventTreatment::discard; } + else if(argn == 5 && std::string(argv[4]) == "split"){ + --argn; + treat[qqxexb] = EventTreatment::reweight; + treat[qqxexf] = EventTreatment::reweight; + treat[FKL] = EventTreatment::discard; + } + else if(argn == 5 && std::string(argv[4]) == "qqxmid"){ + --argn; + treat[qqxmid] = EventTreatment::reweight; + treat[FKL] = EventTreatment::discard; + } if(argn != 4){ std::cerr << "Usage: check_res eventfile xsection tolerance [uno]"; return EXIT_FAILURE; } const double xsec_ref = std::stod(argv[2]); const double tolerance = std::stod(argv[3]); HEJ::istream in{argv[1]}; LHEF::Reader reader{in}; HEJ::PhaseSpacePointConfig psp_conf; psp_conf.jet_param = HEJ::JetParameters{jet_def, jetptmin}; psp_conf.min_extparton_pt = extpartonptmin; psp_conf.max_ext_soft_pt_fraction = max_ext_soft_pt_fraction; HEJ::MatrixElementConfig ME_conf; ME_conf.log_correction = log_corr; ME_conf.Higgs_coupling = HEJ::HiggsCouplingSettings{}; HEJ::EventReweighterConfig conf; conf.psp_config = std::move(psp_conf); conf.ME_config = std::move(ME_conf); conf.jet_param = psp_conf.jet_param; conf.treat = treat; reader.readEvent(); const bool has_Higgs = std::find( begin(reader.hepeup.IDUP), end(reader.hepeup.IDUP), 25 ) != end(reader.hepeup.IDUP); const double mu = has_Higgs?125.:91.188; HEJ::ScaleGenerator scale_gen{ {{std::to_string(mu), HEJ::FixedScale{mu}}}, {}, 1. }; HEJ::Mixmax ran{}; HEJ::EventReweighter hej{reader.heprup, std::move(scale_gen), conf, ran}; double xsec = 0.; double xsec_err = 0.; do{ auto ev_data = HEJ::Event::EventData{reader.hepeup}; ev_data.reconstruct_intermediate(); HEJ::Event ev{ ev_data.cluster( Born_jet_def, Born_jetptmin ) }; auto resummed_events = hej.reweight(ev, 100); for(auto const & ev: resummed_events) { ASSERT(correct_colour(ev)); xsec += ev.central().weight; xsec_err += ev.central().weight*ev.central().weight; } } while(reader.readEvent()); xsec_err = std::sqrt(xsec_err); const double significance = std::abs(xsec - xsec_ref) / std::sqrt( xsec_err*xsec_err + tolerance*tolerance ); std::cout << xsec_ref << " +/- " << tolerance << " ~ " << xsec << " +- " << xsec_err << " => " << significance << " sigma\n"; if(significance > 3.){ std::cerr << "Cross section is off by over 3 sigma!\n"; return EXIT_FAILURE; } }