Page MenuHomeHEPForge
Files F8724466

No OneTemporary
Actions

View Options

diff --git a/t/test_ME_generic.cc b/t/test_ME_generic.cc
index 2168702..02c4556 100644
--- a/t/test_ME_generic.cc
+++ b/t/test_ME_generic.cc
@@ -1,108 +1,120 @@
// Generic tester for the ME for a given set of PSP
// reference weights and PSP (as LHE file) have to be given as _individual_ files
#include <fstream>
+#include <random>
+#include <algorithm>
#include "LHEF/LHEF.h"
#include "HEJ/MatrixElement.hh"
#include "HEJ/Event.hh"
#include "HEJ/YAMLreader.hh"
#include "HEJ/stream.hh"
constexpr double alpha_s = 0.118;
constexpr double ep = 1e-6;
+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);
+}
+
void dump(HEJ::Event const & ev){
{
LHEF::Writer writer{std::cout};
std::cout << std::setprecision(6);
writer.hepeup = to_HEPEUP(std::move(ev), nullptr);
writer.writeEvent();
}
std::cout << "Rapidity ordering:\n";
for(const auto & part: ev.outgoing()){
std::cout << std::setw(2) << part.type << ": "<< std::setw(7) << part.rapidity() << std::endl;
}
}
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};
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()){
++i;
+ HEJ::Event::EventData data{reader.hepeup};
+ shuffle_particles(data);
+
HEJ::Event event{
- HEJ::Event::EventData{reader.hepeup}(
+ data.cluster(
config.resummation_jets.def,
config.resummation_jets.min_pt
)
};
const double our_ME = ME.tree(event).central;
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::cout.precision(16);
std::cout<< "Large difference in PSP " << i << "\nis: "<<our_ME << " should: " << ref_ME << " => difference: " << diff << std::endl;
std::cout.precision(precision);
dump(event);
return EXIT_FAILURE;
}
}
}
wgt_file.close();
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;
}
diff --git a/t/test_classify.cc b/t/test_classify.cc
index aad756b..13081e3 100644
--- a/t/test_classify.cc
+++ b/t/test_classify.cc
@@ -1,51 +1,62 @@
+#include <random>
+#include <algorithm>
+
#include "LHEF/LHEF.h"
#include "HEJ/stream.hh"
#include "HEJ/event_types.hh"
#include "HEJ/Event.hh"
namespace{
constexpr double min_jet_pt = 30.;
const fastjet::JetDefinition jet_def{fastjet::kt_algorithm, 0.4};
using namespace HEJ::event_type;
static const std::vector<EventType> results{
unob,FKL,FKL,FKL,FKL,FKL,FKL,unob,FKL,unob,FKL,FKL,FKL,unof,FKL,unob,FKL,
FKL,unob,unob,FKL,FKL,unob,FKL,FKL,FKL,FKL,FKL,FKL,FKL,FKL,FKL,FKL,FKL,unof,
FKL,FKL,unof,FKL,FKL,FKL,FKL,FKL,unof,FKL,FKL,FKL,unof,FKL,FKL,unob,unof,
FKL,unof,FKL,unob,FKL,FKL,unob,FKL,unob,unof,unob,unof,FKL,FKL,FKL,FKL,FKL,
FKL,FKL,FKL,FKL,FKL,FKL,FKL,unob,FKL,FKL,FKL,FKL,FKL,FKL,FKL,FKL,unob,FKL,
FKL,FKL,FKL,unof,FKL,unob,FKL,FKL,FKL,FKL,FKL,FKL,FKL,FKL,FKL,FKL,unob,FKL,
FKL,FKL,FKL,FKL,unob,FKL,unob,unob,FKL,FKL,FKL,FKL,FKL,FKL,FKL,FKL,FKL,unof,unob,FKL
};
+
+ 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);
+ }
}
int main(int argn, char** argv) {
if(argn != 2){
std::cerr << "Usage: test_classify eventfile";
return EXIT_FAILURE;
}
HEJ::istream in{argv[1]};
LHEF::Reader reader{in};
LHEF::Writer writer{std::cerr};
writer.heprup = reader.heprup;
for(auto const & expected: results){
reader.readEvent();
+ HEJ::Event::EventData data{reader.hepeup};
+ shuffle_particles(data);
const HEJ::Event ev{
- HEJ::Event::EventData{reader.hepeup}(
+ data.cluster(
jet_def, min_jet_pt
)
};
if(ev.type() != expected){
using HEJ::event_type::names;
writer.hepeup = reader.hepeup;
std::cerr << "wrong classification of event:\n";
writer.writeEvent();
std::cerr << "classified as " << names[ev.type()]
<< ", is " << names[expected] << '\n';
return EXIT_FAILURE;
}
}
}
diff --git a/t/test_scale_arithmetics.cc b/t/test_scale_arithmetics.cc
index ca48b3a..b640726 100644
--- a/t/test_scale_arithmetics.cc
+++ b/t/test_scale_arithmetics.cc
@@ -1,88 +1,99 @@
// Generic tester for the ME for a given set of PSP
// reference weights and PSP (as LHE file) have to be given as _individual_ files
#include <fstream>
+#include <random>
+#include <algorithm>
#include "LHEF/LHEF.h"
#include "HEJ/EventReweighter.hh"
#include "HEJ/make_RNG.hh"
#include "HEJ/Event.hh"
#include "HEJ/YAMLreader.hh"
#include "HEJ/stream.hh"
constexpr double alpha_s = 0.118;
constexpr double ep = 1e-13;
void dump(HEJ::Event const & ev){
{
LHEF::Writer writer{std::cout};
std::cout << std::setprecision(6);
writer.hepeup = to_HEPEUP(std::move(ev), nullptr);
writer.writeEvent();
}
std::cout << "Rapidity ordering:\n";
for(const auto & part: ev.outgoing()){
std::cout << std::setw(2) << part.type << ": "<< std::setw(7) << part.rapidity() << std::endl;
}
}
+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);
+}
+
int main(int argn, char** argv){
if(argn != 3){
std::cerr << "\n# Usage:\n."<< argv[0] <<" config.yml input_file.lhe\n\n";
return EXIT_FAILURE;
}
HEJ::Config config = HEJ::load_config(argv[1]);
config.scales = HEJ::to_ScaleConfig(
YAML::Load("scales: [H_T, 1 * H_T, 2/2 * H_T, 2*H_T/2, H_T/2*2, H_T/2/2*4, H_T*H_T/H_T]")
);
HEJ::istream in{argv[2]};
LHEF::Reader reader{in};
auto ran = HEJ::make_RNG(config.rng.name, config.rng.seed);
HEJ::ScaleGenerator scale_gen{
config.scales.base,
config.scales.factors,
config.scales.max_ratio
};
HEJ::EventReweighter resum{
reader.heprup,
std::move(scale_gen),
to_EventReweighterConfig(config),
*ran
};
size_t i = 0;
while(reader.readEvent()){
++i;
+ HEJ::Event::EventData data{reader.hepeup};
+ shuffle_particles(data);
+
HEJ::Event event{
- HEJ::Event::EventData{reader.hepeup}(
+ data.cluster(
config.resummation_jets.def,
config.resummation_jets.min_pt
)
};
auto resummed = resum.reweight(event, config.trials);
for(auto && ev: resummed) {
for(auto &&var: ev.variations()) {
if(std::abs(var.muf - ev.central().muf) > ep) {
std::cerr
<< std::setprecision(15)
<< "unequal scales: " << var.muf
<< " != " << ev.central().muf << '\n'
<< "in resummed event:\n";
dump(ev);
std::cerr << "\noriginal event:\n";
dump(event);
return EXIT_FAILURE;
}
}
}
}
}

File Metadata

Mime Type
text/x-diff
Expires
Mon, Jan 20, 11:32 PM (1 d, 6 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
4242883
Default Alt Text
(8 KB)

Event Timeline