Page MenuHomeHEPForge
Files F8723890

No OneTemporary
Actions

View Options

diff --git a/include/HEJ/EventReader.hh b/include/HEJ/EventReader.hh
new file mode 100644
index 0000000..205e3ff
--- /dev/null
+++ b/include/HEJ/EventReader.hh
@@ -0,0 +1,44 @@
+/** \file
+ * \brief Header file for event reader interface
+ *
+ * This header defines an abstract base class for reading events from files.
+ *
+ * \authors Jeppe Andersen, Marian Heil, Andreas Maier, Jennifer Smillie
+ * \date 2019
+ * \copyright GPLv2 or later
+ */
+#pragma once
+
+#include<string>
+#include<memory>
+
+#include "LHEF/LHEF.h"
+
+namespace HEJ{
+ class EventData;
+
+ // Abstract base class for reading events from files
+ struct EventReader {
+ //! Read an event
+ virtual bool read_event() = 0;
+
+ //! Access header text
+ virtual std::string const & header() const = 0;
+
+ //! Access run information
+ virtual LHEF::HEPRUP const & heprup() const = 0;
+
+ //! Access last read event
+ virtual LHEF::HEPEUP const & hepeup() const = 0;
+
+ virtual ~EventReader() = default;
+ };
+
+ //! Factory function for event readers
+ /**
+ * @param infile The name of the input file
+ * @returns A pointer to an instance of an EventReader
+ * for the input file
+ */
+ std::unique_ptr<EventReader> make_reader(std::string const & filename);
+}
diff --git a/include/HEJ/LesHouchesReader.hh b/include/HEJ/LesHouchesReader.hh
new file mode 100644
index 0000000..e79cc08
--- /dev/null
+++ b/include/HEJ/LesHouchesReader.hh
@@ -0,0 +1,52 @@
+/** \file
+ * \brief Header file for reading events in the Les Houches Event File format.
+ *
+ * \authors Jeppe Andersen, Marian Heil, Andreas Maier, Jennifer Smillie
+ * \date 2019
+ * \copyright GPLv2 or later
+ */
+#pragma once
+
+#include "HEJ/EventReader.hh"
+#include "HEJ/Event.hh"
+#include "HEJ/stream.hh"
+
+#include "LHEF/LHEF.h"
+
+namespace HEJ{
+
+ //! Class for writing events to a file in the Les Houches Event File format
+ class LesHouchesReader : public EventReader{
+ public:
+ //! Contruct object reading from the given file
+ explicit LesHouchesReader(std::string const & filename):
+ stream_{filename},
+ reader_{stream_}
+ {}
+
+ //! Read an event
+ bool read_event() override {
+ return reader_.readEvent();
+ }
+
+ //! Access header text
+ std::string const & header() const override {
+ return reader_.headerBlock;
+ }
+
+ //! Access run information
+ LHEF::HEPRUP const & heprup() const override {
+ return reader_.heprup;
+ }
+
+ //! Access last read event
+ LHEF::HEPEUP const & hepeup() const override {
+ return reader_.hepeup;
+ }
+
+ private:
+ HEJ::istream stream_;
+ LHEF::Reader reader_;
+ };
+
+}
diff --git a/src/EventReader.cc b/src/EventReader.cc
new file mode 100644
index 0000000..ec4c3fd
--- /dev/null
+++ b/src/EventReader.cc
@@ -0,0 +1,9 @@
+#include "HEJ/EventReader.hh"
+#include "HEJ/LesHouchesReader.hh"
+#include "HEJ/utility.hh"
+
+namespace HEJ {
+ std::unique_ptr<EventReader> make_reader(std::string const & filename) {
+ return std::make_unique<LesHouchesReader>(filename);
+ }
+}
diff --git a/src/bin/HEJ.cc b/src/bin/HEJ.cc
index 77d6f1e..d799cc4 100644
--- a/src/bin/HEJ.cc
+++ b/src/bin/HEJ.cc
@@ -1,189 +1,188 @@
/**
* \authors Jeppe Andersen, Tuomas Hapola, Marian Heil, Andreas Maier, Jennifer Smillie
* \date 2019
* \copyright GPLv2 or later
*/
#include <array>
#include <chrono>
#include <iostream>
#include <limits>
#include <memory>
#include <numeric>
#include "yaml-cpp/yaml.h"
-#include "LHEF/LHEF.h"
-
#include "fastjet/ClusterSequence.hh"
#include "HEJ/CombinedEventWriter.hh"
#include "HEJ/config.hh"
#include "HEJ/CrossSectionAccumulator.hh"
#include "HEJ/Event.hh"
+#include "HEJ/EventReader.hh"
#include "HEJ/EventReweighter.hh"
#include "HEJ/get_analysis.hh"
#include "HEJ/make_RNG.hh"
#include "HEJ/ProgressBar.hh"
#include "HEJ/stream.hh"
#include "HEJ/Version.hh"
#include "HEJ/YAMLreader.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));
}
HEJ::Config load_config(char const * filename){
try{
return HEJ::load_config(filename);
}
catch(std::exception const & exc){
std::cerr << "Error: " << exc.what() << '\n';
std::exit(EXIT_FAILURE);
}
}
std::unique_ptr<HEJ::Analysis> get_analysis(
YAML::Node const & parameters
){
try{
return HEJ::get_analysis(parameters);
}
catch(std::exception const & exc){
std::cerr << "Failed to load analysis: " << exc.what() << '\n';
std::exit(EXIT_FAILURE);
}
}
// unique_ptr is a workaround:
// HEJ::optional is a better fit, but gives spurious errors with g++ 7.3.0
std::unique_ptr<HEJ::ProgressBar<double>> make_progress_bar(
std::vector<double> const & xs
) {
if(xs.empty()) return {};
const double Born_xs = std::accumulate(begin(xs), end(xs), 0.);
return std::make_unique<HEJ::ProgressBar<double>>(std::cout, Born_xs);
}
std::string time_to_string(const time_t time){
char s[30];
struct tm * p = localtime(&time);
strftime(s, 30, "%a %b %d %Y %H:%M:%S", p);
return s;
}
int main(int argn, char** argv) {
using clock = std::chrono::system_clock;
if (argn < 3) {
std::cerr << "\n# Usage:\n."<< argv[0] <<" config_file input_file\n\n";
return EXIT_FAILURE;
}
const auto start_time = clock::now();
{
std::cout << "Starting " << HEJ::Version::package_name_full()
<< ", revision " << HEJ::Version::revision() << " ("
<< time_to_string(clock::to_time_t(start_time)) << ")" << std::endl;
}
fastjet::ClusterSequence::print_banner();
// read configuration
const HEJ::Config config = load_config(argv[1]);
-
- HEJ::istream in{argv[2]};
+ auto reader = HEJ::make_reader(argv[2]);
+ assert(reader);
std::unique_ptr<HEJ::Analysis> analysis = get_analysis(
config.analysis_parameters
);
assert(analysis != nullptr);
- LHEF::Reader reader{in};
- auto heprup = reader.heprup;
+ auto heprup = reader->heprup();
heprup.generators.emplace_back(LHEF::XMLTag{});
heprup.generators.back().name = HEJ::Version::package_name();
heprup.generators.back().version = HEJ::Version::String();
HEJ::CombinedEventWriter writer{config.output, std::move(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);
+ const int input_events = event_number(reader->header());
global_reweight = input_events/static_cast<double>(max_events);
std::cout << "Processing " << max_events
<< " out of " << input_events << " events\n";
}
HEJ::ScaleGenerator scale_gen{
config.scales.base,
config.scales.factors,
config.scales.max_ratio
};
auto ran = HEJ::make_RNG(config.rng.name, config.rng.seed);
assert(ran != nullptr);
HEJ::EventReweighter hej{
- reader.heprup,
+ reader->heprup(),
std::move(scale_gen),
to_EventReweighterConfig(config),
*ran
};
int nevent = 0;
std::array<int, HEJ::event_type::last_type + 1>
nevent_type{0}, nfailed_type{0};
- auto progress = make_progress_bar(reader.heprup.XSECUP);
+ auto progress = make_progress_bar(reader->heprup().XSECUP);
HEJ::CrossSectionAccumulator xs;
// Loop over the events in the inputfile
- while(reader.readEvent()){
+ while(reader->read_event()){
// reweight events so that the total cross section is conserved
- reader.hepeup.setWeight(0, global_reweight * reader.hepeup.weight());
+ auto hepeup = reader->hepeup();
+ hepeup.setWeight(0, global_reweight * hepeup.weight());
if(nevent == max_events) break;
++nevent;
// calculate HEJ weight
HEJ::Event FO_event{
- HEJ::Event::EventData{reader.hepeup}(
+ HEJ::Event::EventData{hepeup}(
config.fixed_order_jets.def, config.fixed_order_jets.min_pt
)
};
auto resummed_events = hej.reweight(FO_event, config.trials);
++nevent_type[FO_event.type()];
if(resummed_events.empty()) ++nfailed_type[FO_event.type()];
for(auto const & ev: resummed_events){
//TODO: move pass_cuts to after phase space point generation
if(analysis->pass_cuts(ev, FO_event)){
analysis->fill(ev, FO_event);
writer.write(ev);
xs.fill(ev);
}
}
if(progress) progress->increment(FO_event.central().weight);
} // main event loop
std::cout << '\n';
analysis->finalise();
using namespace HEJ::event_type;
std::cout<< "Events processed: " << nevent << '\n';
for(size_t ev_type = first_type; ev_type <= last_type; ++ev_type){
std::cout << '\t' << names[ev_type] << ": " << nevent_type[ev_type]
<< ", failed to reconstruct " << nfailed_type[ev_type]
<< '\n';
}
std::cout << '\n' << xs << '\n';
std::chrono::duration<double> run_time = (clock::now() - start_time);
std::cout << "Finished " << HEJ::Version::package_name() << " at "
<< time_to_string(clock::to_time_t(clock::now()))
<< "\n=> Runtime: " << run_time.count() << " sec ("
<< nevent/run_time.count() << " Events/sec).\n";
}

File Metadata

Mime Type
text/x-diff
Expires
Mon, Jan 20, 9:50 PM (1 d, 4 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
4242539
Default Alt Text
(9 KB)

Event Timeline