Page MenuHomeHEPForge
Files F9501390

No OneTemporary
Actions

View Options

diff --git a/t/test_ME_generic.cc b/t/test_ME_generic.cc
index 342fe79..290fc1c 100644
--- a/t/test_ME_generic.cc
+++ b/t/test_ME_generic.cc
@@ -1,134 +1,132 @@
/**
* \brief Generic tester for the ME for a given set of PSP
*
* \note reference weights and PSP (as LHE file) have to be given as
* _individual_ files
*
* \authors The HEJ collaboration (see AUTHORS for details)
* \date 2019
* \copyright GPLv2 or later
*/
#include <algorithm>
#include <fstream>
#include <math.h>
#include <random>
-#include "LHEF/LHEF.h"
-
#include "HEJ/Event.hh"
+#include "HEJ/EventReader.hh"
#include "HEJ/MatrixElement.hh"
#include "HEJ/stream.hh"
#include "HEJ/YAMLreader.hh"
constexpr double alpha_s = 0.118;
constexpr double ep = 1e-5;
void shuffle_particles(HEJ::Event::EventData & ev) {
static std::mt19937_64 ran{0};
std::shuffle(begin(ev.incoming), end(ev.incoming), ran);
std::shuffle(begin(ev.outgoing), end(ev.outgoing), ran);
}
enum MEComponent {tree, virt};
MEComponent guess_component(std::string const & data_file) {
if(data_file.find("virt") != data_file.npos) return MEComponent::virt;
return MEComponent::tree;
}
int main(int argn, char** argv){
if(argn != 4 && argn != 5){
std::cerr << "\n# Usage:\n."<< argv[0] <<" config.yml ME_weights input_file.lhe\n\n";
return EXIT_FAILURE;
}
bool OUTPUT_MODE = false;
if(argn == 5 && std::string("OUTPUT")==std::string(argv[4]))
OUTPUT_MODE = true;
const HEJ::Config config = HEJ::load_config(argv[1]);
std::fstream wgt_file;
if ( OUTPUT_MODE ) {
std::cout << "_______________________USING OUTPUT MODE!_______________________" << std::endl;
wgt_file.open(argv[2], std::fstream::out);
wgt_file.precision(10);
} else {
wgt_file.open(argv[2], std::fstream::in);
}
- HEJ::istream in{argv[3]};
- LHEF::Reader reader{in};
+ auto reader{ HEJ::make_reader(argv[3])};
const auto component = guess_component(argv[2]);
HEJ::MatrixElement ME{
[](double){ return alpha_s; },
HEJ::to_MatrixElementConfig(config)
};
double max_ratio = 0.;
size_t idx_max_ratio = 0;
HEJ::Event ev_max_ratio(HEJ::Event::EventData{}.cluster(
config.resummation_jets.def,0
)
);
double av_ratio = 0;
size_t i = 0;
- while(reader.readEvent()){
+ while(reader->read_event()){
++i;
- HEJ::Event::EventData data{reader.hepeup};
+ HEJ::Event::EventData data{reader->hepeup()};
shuffle_particles(data);
HEJ::Event event{
data.cluster(
config.resummation_jets.def,
config.resummation_jets.min_pt
)
};
const double our_ME = (component == MEComponent::tree)?
ME.tree(event).central:
ME.virtual_corrections(event).central
;
if(!isfinite(our_ME)){
std::cerr << "Found non-finite ME ("<< our_ME <<")\n" << event << std::endl;
return EXIT_FAILURE;
}
if ( OUTPUT_MODE ) {
wgt_file << our_ME << std::endl;
} else {
std::string line;
if(!std::getline(wgt_file,line)) break;
const double ref_ME = std::stod(line);
const double diff = std::abs(our_ME/ref_ME-1.);
av_ratio+=diff;
if( diff > max_ratio ) {
max_ratio = diff;
idx_max_ratio = i;
ev_max_ratio = event;
}
if( diff > ep ){
size_t precision(std::cout.precision());
std::cerr.precision(16);
std::cerr<< "Large difference in PSP " << i << "\nis: "<<our_ME
<< " should: " << ref_ME << " => difference: " << diff
<< event << std::endl;
std::cerr.precision(precision);
return EXIT_FAILURE;
}
}
}
wgt_file.close();
if ( i<100 )
throw std::invalid_argument{"Not enough PSP tested"};
if ( !OUTPUT_MODE ) {
size_t precision(std::cout.precision());
std::cout.precision(16);
std::cout << "Avg ratio after " << i << " PSP: " << av_ratio/i << std::endl;
std::cout << "maximal ratio at " << idx_max_ratio << ": " << max_ratio << std::endl;
std::cout.precision(precision);
}
return EXIT_SUCCESS;
}

File Metadata

Mime Type
text/x-diff
Expires
Sun, Feb 23, 2:25 PM (15 h, 34 m)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
4482603
Default Alt Text
(4 KB)

Event Timeline