From edd991d9aee6613a19033ff9c08f68ab07e71b29 Mon Sep 17 00:00:00 2001 From: Marco Giacalone Date: Sun, 12 Oct 2025 12:38:35 +0200 Subject: [PATCH 1/2] Added first working version of HERWIG7 external generator --- MC/config/examples/herwig/LHC.in | 77 ++++ MC/config/examples/herwig/generator_Herwig.C | 405 +++++++++++++++++++ MC/config/examples/ini/GeneratorHERWIG7.ini | 3 + 3 files changed, 485 insertions(+) create mode 100644 MC/config/examples/herwig/LHC.in create mode 100644 MC/config/examples/herwig/generator_Herwig.C create mode 100644 MC/config/examples/ini/GeneratorHERWIG7.ini diff --git a/MC/config/examples/herwig/LHC.in b/MC/config/examples/herwig/LHC.in new file mode 100644 index 000000000..c189cd8ed --- /dev/null +++ b/MC/config/examples/herwig/LHC.in @@ -0,0 +1,77 @@ +# -*- ThePEG-repository -*- + +################################################################################ +# This file contains our best tune to UE data from ATLAS at 7 TeV. More recent +# tunes and tunes for other centre-of-mass energies as well as more usage +# instructions can be obtained from this Herwig wiki page: +# http://projects.hepforge.org/herwig/trac/wiki/MB_UE_tunes +# The model for soft interactions and diffractions is explained in +# [S. Gieseke, P. Kirchgaesser, F. Loshaj, arXiv:1612.04701] +################################################################################ + +read snippets/PPCollider.in + +################################################## +# Technical parameters for this run +################################################## +cd /Herwig/Generators +################################################## +# LHC physics parameters (override defaults here) +################################################## +set EventGenerator:EventHandler:LuminosityFunction:Energy 13000.0 + + +# Minimum Bias +read snippets/MB.in + +# Recommended set of parameters for MB/UE runs + +set /Herwig/Hadronization/ColourReconnector:ReconnectionProbability 0.5 +set /Herwig/UnderlyingEvent/MPIHandler:pTmin0 3.502683 +set /Herwig/UnderlyingEvent/MPIHandler:InvRadius 1.402055 +set /Herwig/UnderlyingEvent/MPIHandler:Power 0.416852 +set /Herwig/Partons/RemnantDecayer:ladderPower -0.08 +set /Herwig/Partons/RemnantDecayer:ladderNorm 0.95 + +# Diffraction model +# read snippets/Diffraction.in + +# Read in snippet in order to use baryonic reconnection model with modified gluon splitting (uds) +# For more details see [S. Gieseke, P. Kirchgaeßer, S. Plätzer. arXiv:1710.10906]] +############################################################################################## + +# read BaryonicReconnection.in + + +################################################## +# Analyses +################################################## + + +# cd /Herwig/Analysis +# create ThePEG::RivetAnalysis RivetAnalysis RivetAnalysis.so + +# cd /Herwig/Generators +# insert EventGenerator:AnalysisHandlers 0 /Herwig/Analysis/RivetAnalysis + +# insert /Herwig/Analysis/RivetAnalysis:Analyses 0 ATLAS_20XX_XXXXXXX + + +#set /Herwig/Analysis/Plot:EventNumber 54 +#cd /Herwig/Generators +#insert EventGenerator:AnalysisHandlers 0 /Herwig/Analysis/Plot + +# read /Herwig/snippets/HepMC.in +#create ThePEG::HepMCFile /Herwig/Analysis/HepMC HepMCAnalysis.so +#set /Herwig/Analysis/HepMC:PrintEvent $NEV +#set /Herwig/Analysis/HepMC:Format GenEvent +#set /Herwig/Analysis/HepMC:Units GeV_mm +#set /Herwig/Analysis/HepMC:Filename ${FIFOPATH} +#insert /Herwig/Generators/EventGenerator:AnalysisHandlers 0 /Herwig/Analysis/HepMC + + +################################################## +# Save run for later usage with 'Herwig run' +################################################## +# cd /Herwig/Generators +saverun LHC EventGenerator \ No newline at end of file diff --git a/MC/config/examples/herwig/generator_Herwig.C b/MC/config/examples/herwig/generator_Herwig.C new file mode 100644 index 000000000..b5e540f8d --- /dev/null +++ b/MC/config/examples/herwig/generator_Herwig.C @@ -0,0 +1,405 @@ +// Set include paths for ThePEG headers +R__ADD_INCLUDE_PATH($THEPEG_ROOT/include/) +R__ADD_INCLUDE_PATH($THEPEG_ROOT/../../GSL/latest/include/) +R__ADD_INCLUDE_PATH($HEPMC3_ROOT/include/) +R__ADD_INCLUDE_PATH($HERWIG_ROOT/include/) +#define SKIP_HEPMC_CONVERSION 1 +#define HAVE_HEPMC3 1 +#define ThePEG_DEBUG_LEVEL 1 + +// O2DPG and ROOT includes +#include "FairGenerator.h" +#include "FairPrimaryGenerator.h" +#include "fairlogger/Logger.h" +#include "TRandom3.h" +#include "TParticle.h" +#include "TParticlePDG.h" +#include "TDatabasePDG.h" +#include "TLorentzVector.h" +#include "TMath.h" +#include +#include +#include +#include +#include + +// Undefine conflicting macros before including ThePEG headers +#ifdef B0 +#undef B0 +#endif +#ifdef B50 +#undef B50 +#endif +#ifdef B75 +#undef B75 +#endif +#ifdef B110 +#undef B110 +#endif +#ifdef B134 +#undef B134 +#endif +#ifdef B150 +#undef B150 +#endif +#ifdef B200 +#undef B200 +#endif +#ifdef B300 +#undef B300 +#endif +#ifdef B600 +#undef B600 +#endif +#ifdef B1200 +#undef B1200 +#endif +#ifdef B1800 +#undef B1800 +#endif +#ifdef B2400 +#undef B2400 +#endif +#ifdef B4800 +#undef B4800 +#endif +#ifdef B9600 +#undef B9600 +#endif +#ifdef B19200 +#undef B19200 +#endif +#ifdef B38400 +#undef B38400 +#endif + +// ThePEG and Herwig includes +#include "ThePEG/Repository/EventGenerator.h" +#include "ThePEG/EventRecord/Event.h" +#include "ThePEG/EventRecord/Particle.h" +#include "ThePEG/EventRecord/Step.h" +#include "ThePEG/Config/ThePEG.h" +#include "ThePEG/PDT/ParticleData.h" +#include "ThePEG/Vectors/HepMCConverter.h" +#include "ThePEG/Repository/Repository.h" +#include "ThePEG/Repository/BaseRepository.h" +#include "ThePEG/Utilities/DynamicLoader.h" +#include "ThePEG/Persistency/PersistentIStream.h" + +// Herwig specific includes +#include "Herwig/API/HerwigAPI.h" +#include "Herwig/API/HerwigUI.h" + +// Subclass of HerwigUI to provide minimal implementation for our use case +class SimpleHerwigUI : public Herwig::HerwigUI +{ +public: + SimpleHerwigUI(const std::string &inFile, + Herwig::RunMode::Mode mode = Herwig::RunMode::READ) + : m_inFile(inFile), m_mode(mode), + m_in(inFile), m_out(std::cout), m_err(std::cerr) + { + if (!m_in) + throw std::runtime_error("Cannot open Herwig input file: " + inFile); + Dirs.reserve(5); + std::string hDir = std::getenv("HERWIG_ROOT"); + if (!hDir.empty()) + Dirs.push_back(hDir + "/share/Herwig"); + } + + Herwig::RunMode::Mode runMode() const override { return m_mode; } + + std::string repository() const override { + std::string rpo = std::getenv("HERWIG_ROOT"); + rpo.append("/share/Herwig/HerwigDefaults.rpo"); + return rpo; + } + std::string inputfile() const override { return m_inFile; } + std::string setupfile() const override { return ""; } + + bool resume() const override { return false; } + bool tics() const override { return false; } + std::string tag() const override { return ""; } + std::string integrationList() const override { return ""; } + + const std::vector &prependReadDirectories() const override + { + return Dirs; + } + const std::vector &appendReadDirectories() const override + { + static std::vector empty; + return empty; + } + + long N() const override { return 10; } // number of events + int seed() const override { return 1234; } + int jobs() const override { return 1; } + unsigned int jobSize() const override { return 1; } + unsigned int maxJobs() const override { return 1; } + + void quitWithHelp() const override { std::exit(1); } + void quit() const override { std::exit(1); } + + std::ostream &outStream() const override { return m_out; } + std::ostream &errStream() const override { return m_err; } + std::istream &inStream() const override { return m_in; } + +private: + std::string m_inFile; + Herwig::RunMode::Mode m_mode; + mutable std::ifstream m_in; + std::ostream &m_out; + std::ostream &m_err; + std::vector Dirs; +}; + +namespace o2 +{ +namespace eventgen +{ + +/// HERWIG7 event generator using ThePEG interface +/// Author: Marco Giacalone (marco.giacalone@cern.ch) +/// Based on the O2DPG external generator patterns +class GeneratorHerwig : public Generator +{ +public: + /// Default constructor + GeneratorHerwig(const std::string& configFile = "LHC.in") + : mConfigFile(configFile) + , mEventGenerator(nullptr) + { + LOG(info) << "HERWIG7 Generator initialized"; + LOG(info) << "Config file: " << mConfigFile; + std::string extension = mConfigFile.substr(mConfigFile.find_last_of(".")); + if( extension == ".in" ) { + mIsInputFile = true; + LOG(info) << "Using input file for configuration"; + } else if(std::find(mConfigFile.end()-4, mConfigFile.end(), '.run') != mConfigFile.end()) { + mIsInputFile = false; + LOG(info) << "Using run file for configuration"; + } else { + LOG(fatal) << "No file extension found in config file: " << mConfigFile; + exit(1); + } + } + + /// Destructor + ~GeneratorHerwig(){ + delete mEventGenerator; + }; + + /// Initialize the generator + Bool_t Init() override + { + LOG(info) << "Initializing HERWIG7 Generator"; + + try { + if (mIsInputFile) { + // Process .in file + return initFromInputFile(); + } else { + // Process .run file directly + return initFromRunFile(); + } + + } catch (const std::exception& e) { + LOG(fatal) << "Exception during HERWIG7 initialization: " << e.what(); + return kFALSE; + } + } + + /// Generate event + Bool_t generateEvent() override + { + if (!mEventGenerator) { + LOG(error) << "Event generator not initialized"; + return kFALSE; + } + try { + // Clear previous event particles + hParticles.clear(); + + // Generate the event + ThePEG::EventPtr event = mEventGenerator->shoot(); + + if (!event) { + LOG(error) << "Failed to generate event"; + return kFALSE; + } + + // Convert ThePEG event to TParticle format + convertEvent(event); + LOG(debug) << "Herwig7 generated " << hParticles.size(); + + return kTRUE; + + } catch (const std::exception& e) { + LOG(error) << "Exception during event generation: " << e.what(); + return kFALSE; + } + } + + /// Import particles for transport + Bool_t importParticles() override + { + if (hParticles.empty()) { + LOG(warning) << "No particles to import"; + return kFALSE; + } + + // Add particles to the primary generator + for (const auto& particle : hParticles) { + mParticles.push_back(particle); + } + + return kTRUE; + } + +private: + std::string mConfigFile; ///< HERWIG config file (.in or .run) + Bool_t mIsInputFile = false; ///< True for .in files, false for .run files + ThePEG::EGPtr mEventGenerator; ///< ThePEG event generator + std::vector hParticles; ///< Generated Herwig particles + + void printHerwigSearchPaths() + { + const auto &list = ThePEG::Repository::listReadDirs(); + + LOG(info) << "Append directories:\n"; + for (const auto &p : list) + LOG(info) << " " << p << "\n"; + } + + /// Initialize from .in file + Bool_t initFromInputFile() + { + LOG(info) << "Initializing from .in file: " << mConfigFile; + + try { + using namespace ThePEG; + SimpleHerwigUI ui(mConfigFile, Herwig::RunMode::READ); + + Herwig::API::read(ui); + printHerwigSearchPaths(); + std::string runFile = mConfigFile; + size_t pos = runFile.find_last_of('.'); + runFile.replace(pos, 4, ".run"); + mConfigFile = runFile; + LOG(info) << "Generated run file: " << runFile; + auto res = initFromRunFile(); + if (!res) { + LOG(error) << "Failed to initialize from generated run file"; + return kFALSE; + } + return kTRUE; + } + catch (const std::exception &e) + { + std::cerr << "Exception: " << e.what() << std::endl; + return kFALSE; + } + } + + /// Initialize from .run file + Bool_t initFromRunFile() + { + LOG(info) << "Initializing from .run file: " << mConfigFile; + + try { + using namespace ThePEG; + + if (!std::ifstream(mConfigFile)) + { + LOG(info) << "Run file does not exist: " << mConfigFile; + return kFALSE; + } + SimpleHerwigUI runui(mConfigFile, Herwig::RunMode::RUN); + // Prepare the generator + mEventGenerator = Herwig::API::prepareRun(runui); + if (!mEventGenerator) + { + std::cerr << "Error: prepareRun() returned null." << std::endl; + return kFALSE; + } + + mEventGenerator->initialize(); + std::cout << "Herwig generator initialized successfully." << std::endl; + return kTRUE; + } + catch (const std::exception &ex) { + LOG(error) << "Standard exception: " << ex.what(); + return kFALSE; + } + } + + /// Convert ThePEG event to TParticle format + void convertEvent(ThePEG::EventPtr event) + { + if (!event) return; + + // Get all particles from the event + const ThePEG::tPVector& particles = event->getFinalState(); + + for (size_t i = 0; i < particles.size(); ++i) { + ThePEG::tPPtr particle = particles[i]; + if (!particle) continue; + + // Get particle properties + int pdgCode = particle->id(); + int status = getFinalStateStatus(particle); + + // Get 4-momentum + ThePEG::LorentzMomentum momentum = particle->momentum(); + double px = momentum.x() / ThePEG::GeV; // Convert to GeV + double py = momentum.y() / ThePEG::GeV; + double pz = momentum.z() / ThePEG::GeV; + double e = momentum.e() / ThePEG::GeV; + + // Get production vertex + const ThePEG::LorentzPoint &vertex = particle->vertex(); + double vx = vertex.x() / ThePEG::mm; // Convert to mm + double vy = vertex.y() / ThePEG::mm; + double vz = vertex.z() / ThePEG::mm; + double vt = vertex.t() / ThePEG::mm; // Convert to mm/c + + // Create TParticle + TParticle tparticle( + pdgCode, status, + -1, -1, -1, -1, // mother and daughter indices (to be set properly) + px, py, pz, e, + vx, vy, vz, vt + ); + tparticle.SetStatusCode(o2::mcgenstatus::MCGenStatusEncoding(tparticle.GetStatusCode(), 0).fullEncoding); + tparticle.SetBit(ParticleStatus::kToBeDone, // + o2::mcgenstatus::getHepMCStatusCode(tparticle.GetStatusCode()) == 1); + + hParticles.push_back(tparticle); + } + + LOG(debug) << "Converted " << hParticles.size() << " particles from HERWIG7 event"; + } + + /// Determine final state status for particle + int getFinalStateStatus(ThePEG::tPPtr particle) + { + // In HERWIG/ThePEG, check if particle is stable + if (particle->children().empty()) { + return 1; // Final state particle + } else { + return 2; // Intermediate particle + } + } +}; + +} // namespace eventgen +} // namespace o2 + +/// Factory function to create HERWIG7 generator from .in file +/// @param inputFile HERWIG input/run file (e.g., "LHC.in"/"LHC.run") +FairGenerator* generateHerwig7(const std::string& inputFile = "LHC.in") +{ + auto generator = new o2::eventgen::GeneratorHerwig(inputFile); + return generator; +} \ No newline at end of file diff --git a/MC/config/examples/ini/GeneratorHERWIG7.ini b/MC/config/examples/ini/GeneratorHERWIG7.ini new file mode 100644 index 000000000..1d31652fc --- /dev/null +++ b/MC/config/examples/ini/GeneratorHERWIG7.ini @@ -0,0 +1,3 @@ +[GeneratorExternal] +fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/examples/epos4herwig/generator_Herwig.C +funcName=generateHerwig7() \ No newline at end of file From 6c68e93d0783a43ebfdf4f7a1fb1f22cb41fe45f Mon Sep 17 00:00:00 2001 From: Marco Giacalone Date: Tue, 14 Oct 2025 17:11:42 +0200 Subject: [PATCH 2/2] Cleaned and improved generator Header files folders are now directly included in the ROOT_INCLUDE_PATHs via the recipes --- MC/config/examples/herwig/LHC.in | 43 +-- MC/config/examples/herwig/generator_Herwig.C | 254 +++++++----------- MC/config/examples/ini/GeneratorHERWIG7.ini | 4 +- .../examples/ini/tests/GeneratorHERWIG7.C | 38 +++ 4 files changed, 132 insertions(+), 207 deletions(-) create mode 100644 MC/config/examples/ini/tests/GeneratorHERWIG7.C diff --git a/MC/config/examples/herwig/LHC.in b/MC/config/examples/herwig/LHC.in index c189cd8ed..cb1ecb4fb 100644 --- a/MC/config/examples/herwig/LHC.in +++ b/MC/config/examples/herwig/LHC.in @@ -18,14 +18,11 @@ cd /Herwig/Generators ################################################## # LHC physics parameters (override defaults here) ################################################## -set EventGenerator:EventHandler:LuminosityFunction:Energy 13000.0 - - +set EventGenerator:EventHandler:LuminosityFunction:Energy 13600.0 # Minimum Bias read snippets/MB.in # Recommended set of parameters for MB/UE runs - set /Herwig/Hadronization/ColourReconnector:ReconnectionProbability 0.5 set /Herwig/UnderlyingEvent/MPIHandler:pTmin0 3.502683 set /Herwig/UnderlyingEvent/MPIHandler:InvRadius 1.402055 @@ -33,45 +30,7 @@ set /Herwig/UnderlyingEvent/MPIHandler:Power 0.416852 set /Herwig/Partons/RemnantDecayer:ladderPower -0.08 set /Herwig/Partons/RemnantDecayer:ladderNorm 0.95 -# Diffraction model -# read snippets/Diffraction.in - -# Read in snippet in order to use baryonic reconnection model with modified gluon splitting (uds) -# For more details see [S. Gieseke, P. Kirchgaeßer, S. Plätzer. arXiv:1710.10906]] -############################################################################################## - -# read BaryonicReconnection.in - - -################################################## -# Analyses -################################################## - - -# cd /Herwig/Analysis -# create ThePEG::RivetAnalysis RivetAnalysis RivetAnalysis.so - -# cd /Herwig/Generators -# insert EventGenerator:AnalysisHandlers 0 /Herwig/Analysis/RivetAnalysis - -# insert /Herwig/Analysis/RivetAnalysis:Analyses 0 ATLAS_20XX_XXXXXXX - - -#set /Herwig/Analysis/Plot:EventNumber 54 -#cd /Herwig/Generators -#insert EventGenerator:AnalysisHandlers 0 /Herwig/Analysis/Plot - -# read /Herwig/snippets/HepMC.in -#create ThePEG::HepMCFile /Herwig/Analysis/HepMC HepMCAnalysis.so -#set /Herwig/Analysis/HepMC:PrintEvent $NEV -#set /Herwig/Analysis/HepMC:Format GenEvent -#set /Herwig/Analysis/HepMC:Units GeV_mm -#set /Herwig/Analysis/HepMC:Filename ${FIFOPATH} -#insert /Herwig/Generators/EventGenerator:AnalysisHandlers 0 /Herwig/Analysis/HepMC - - ################################################## # Save run for later usage with 'Herwig run' ################################################## -# cd /Herwig/Generators saverun LHC EventGenerator \ No newline at end of file diff --git a/MC/config/examples/herwig/generator_Herwig.C b/MC/config/examples/herwig/generator_Herwig.C index b5e540f8d..1c3da112c 100644 --- a/MC/config/examples/herwig/generator_Herwig.C +++ b/MC/config/examples/herwig/generator_Herwig.C @@ -1,110 +1,48 @@ -// Set include paths for ThePEG headers -R__ADD_INCLUDE_PATH($THEPEG_ROOT/include/) -R__ADD_INCLUDE_PATH($THEPEG_ROOT/../../GSL/latest/include/) -R__ADD_INCLUDE_PATH($HEPMC3_ROOT/include/) -R__ADD_INCLUDE_PATH($HERWIG_ROOT/include/) #define SKIP_HEPMC_CONVERSION 1 #define HAVE_HEPMC3 1 -#define ThePEG_DEBUG_LEVEL 1 // O2DPG and ROOT includes #include "FairGenerator.h" #include "FairPrimaryGenerator.h" #include "fairlogger/Logger.h" -#include "TRandom3.h" #include "TParticle.h" -#include "TParticlePDG.h" -#include "TDatabasePDG.h" #include "TLorentzVector.h" -#include "TMath.h" -#include #include #include #include -#include - -// Undefine conflicting macros before including ThePEG headers -#ifdef B0 -#undef B0 -#endif -#ifdef B50 -#undef B50 -#endif -#ifdef B75 -#undef B75 -#endif -#ifdef B110 -#undef B110 -#endif -#ifdef B134 -#undef B134 -#endif -#ifdef B150 -#undef B150 -#endif -#ifdef B200 -#undef B200 -#endif -#ifdef B300 -#undef B300 -#endif -#ifdef B600 -#undef B600 -#endif -#ifdef B1200 -#undef B1200 -#endif -#ifdef B1800 -#undef B1800 -#endif -#ifdef B2400 -#undef B2400 -#endif -#ifdef B4800 -#undef B4800 -#endif -#ifdef B9600 -#undef B9600 -#endif -#ifdef B19200 -#undef B19200 -#endif -#ifdef B38400 -#undef B38400 -#endif - -// ThePEG and Herwig includes + +// ThePEG #include "ThePEG/Repository/EventGenerator.h" #include "ThePEG/EventRecord/Event.h" #include "ThePEG/EventRecord/Particle.h" #include "ThePEG/EventRecord/Step.h" #include "ThePEG/Config/ThePEG.h" #include "ThePEG/PDT/ParticleData.h" -#include "ThePEG/Vectors/HepMCConverter.h" #include "ThePEG/Repository/Repository.h" -#include "ThePEG/Repository/BaseRepository.h" -#include "ThePEG/Utilities/DynamicLoader.h" -#include "ThePEG/Persistency/PersistentIStream.h" -// Herwig specific includes +// Herwig #include "Herwig/API/HerwigAPI.h" #include "Herwig/API/HerwigUI.h" -// Subclass of HerwigUI to provide minimal implementation for our use case +// Subclass of HerwigUI to provide minimal implementation of the abstract class class SimpleHerwigUI : public Herwig::HerwigUI { public: SimpleHerwigUI(const std::string &inFile, - Herwig::RunMode::Mode mode = Herwig::RunMode::READ) + Herwig::RunMode::Mode mode = Herwig::RunMode::READ, int seed = 0) : m_inFile(inFile), m_mode(mode), - m_in(inFile), m_out(std::cout), m_err(std::cerr) + m_in(inFile), m_out(std::cout), m_err(std::cerr), mSeed(seed) { if (!m_in) - throw std::runtime_error("Cannot open Herwig input file: " + inFile); - Dirs.reserve(5); + { + LOG(fatal) << "Cannot open Herwig input file: " << inFile; + exit(1); + } std::string hDir = std::getenv("HERWIG_ROOT"); if (!hDir.empty()) + { Dirs.push_back(hDir + "/share/Herwig"); + } } Herwig::RunMode::Mode runMode() const override { return m_mode; } @@ -132,8 +70,8 @@ public: return empty; } - long N() const override { return 10; } // number of events - int seed() const override { return 1234; } + long N() const override { return 1; } // number of events + int seed() const override { return mSeed; } int jobs() const override { return 1; } unsigned int jobSize() const override { return 1; } unsigned int maxJobs() const override { return 1; } @@ -152,6 +90,7 @@ private: std::ostream &m_out; std::ostream &m_err; std::vector Dirs; + int mSeed = 0; }; namespace o2 @@ -161,52 +100,51 @@ namespace eventgen /// HERWIG7 event generator using ThePEG interface /// Author: Marco Giacalone (marco.giacalone@cern.ch) -/// Based on the O2DPG external generator patterns +/// Based on the O2DPG external generator configurations class GeneratorHerwig : public Generator { public: /// Default constructor - GeneratorHerwig(const std::string& configFile = "LHC.in") + GeneratorHerwig(const std::string& configFile = "LHC.in", int seed = -1) : mConfigFile(configFile) , mEventGenerator(nullptr) { - LOG(info) << "HERWIG7 Generator initialized"; + LOG(info) << "HERWIG7 Generator construction"; LOG(info) << "Config file: " << mConfigFile; std::string extension = mConfigFile.substr(mConfigFile.find_last_of(".")); if( extension == ".in" ) { mIsInputFile = true; LOG(info) << "Using input file for configuration"; - } else if(std::find(mConfigFile.end()-4, mConfigFile.end(), '.run') != mConfigFile.end()) { + } else if(extension == ".run") { mIsInputFile = false; LOG(info) << "Using run file for configuration"; } else { LOG(fatal) << "No file extension found in config file: " << mConfigFile; exit(1); } + if (seed < 0) + { + auto &conf = o2::conf::SimConfig::Instance(); + mSeed = gRandom->Integer(conf.getStartSeed()); + } else { + mSeed = seed; + } + LOG(info) << "Using seed: " << mSeed << " for HERWIG simulation"; } /// Destructor - ~GeneratorHerwig(){ - delete mEventGenerator; - }; + ~GeneratorHerwig() = default; /// Initialize the generator Bool_t Init() override { LOG(info) << "Initializing HERWIG7 Generator"; - - try { - if (mIsInputFile) { - // Process .in file - return initFromInputFile(); - } else { - // Process .run file directly - return initFromRunFile(); - } - - } catch (const std::exception& e) { - LOG(fatal) << "Exception during HERWIG7 initialization: " << e.what(); - return kFALSE; + if (mIsInputFile) { + // Process .in file + return initFromInputFile(); + } else { + // Process .run file directly + return initFromRunFile(); } } @@ -217,28 +155,23 @@ public: LOG(error) << "Event generator not initialized"; return kFALSE; } - try { - // Clear previous event particles - hParticles.clear(); - - // Generate the event - ThePEG::EventPtr event = mEventGenerator->shoot(); - - if (!event) { - LOG(error) << "Failed to generate event"; - return kFALSE; - } - - // Convert ThePEG event to TParticle format - convertEvent(event); - LOG(debug) << "Herwig7 generated " << hParticles.size(); + // Clear previous event particles + hParticles.clear(); - return kTRUE; - - } catch (const std::exception& e) { - LOG(error) << "Exception during event generation: " << e.what(); + // Generate the event + mPEGEvent = mEventGenerator->shoot(); + + if (!mPEGEvent) + { + LOG(error) << "Failed to generate event"; return kFALSE; } + + // Convert ThePEG event to TParticle format + convertEvent(mPEGEvent); + LOG(debug) << "Herwig7 generated " << hParticles.size() << " particles"; + + return kTRUE; } /// Import particles for transport @@ -259,9 +192,11 @@ public: private: std::string mConfigFile; ///< HERWIG config file (.in or .run) - Bool_t mIsInputFile = false; ///< True for .in files, false for .run files + Bool_t mIsInputFile = false; ///< True for .in files, false for .run files ThePEG::EGPtr mEventGenerator; ///< ThePEG event generator std::vector hParticles; ///< Generated Herwig particles + ThePEG::EventPtr mPEGEvent; ///< Pointer to Current event + int mSeed = 0; ///< Random seed for Herwig void printHerwigSearchPaths() { @@ -277,29 +212,27 @@ private: { LOG(info) << "Initializing from .in file: " << mConfigFile; - try { - using namespace ThePEG; - SimpleHerwigUI ui(mConfigFile, Herwig::RunMode::READ); - - Herwig::API::read(ui); - printHerwigSearchPaths(); - std::string runFile = mConfigFile; - size_t pos = runFile.find_last_of('.'); - runFile.replace(pos, 4, ".run"); - mConfigFile = runFile; - LOG(info) << "Generated run file: " << runFile; - auto res = initFromRunFile(); - if (!res) { - LOG(error) << "Failed to initialize from generated run file"; - return kFALSE; - } - return kTRUE; - } - catch (const std::exception &e) - { - std::cerr << "Exception: " << e.what() << std::endl; + using namespace ThePEG; + SimpleHerwigUI ui(mConfigFile, Herwig::RunMode::READ, mSeed); + Herwig::API::read(ui); + // For debugging, print the search paths + // printHerwigSearchPaths(); + // Currently the .run filename is set inside the .in file itself with + // the line "saverun LHC EventGenerator" or similar. We assume this is the same as + // the .in file name with .run extension, so change that string accordingly in your .in files. + std::string runFile = mConfigFile; + size_t pos = runFile.find_last_of('.'); + runFile.replace(pos, 4, ".run"); + pos = runFile.find_last_of('/'); + runFile = (pos != std::string::npos) ? runFile.substr(pos + 1) : runFile; + mConfigFile = runFile; + LOG(info) << "Generated run file: " << runFile; + auto res = initFromRunFile(); + if (!res) { + LOG(error) << "Failed to initialize from generated run file"; return kFALSE; } + return kTRUE; } /// Initialize from .run file @@ -307,31 +240,25 @@ private: { LOG(info) << "Initializing from .run file: " << mConfigFile; - try { - using namespace ThePEG; - - if (!std::ifstream(mConfigFile)) - { - LOG(info) << "Run file does not exist: " << mConfigFile; - return kFALSE; - } - SimpleHerwigUI runui(mConfigFile, Herwig::RunMode::RUN); - // Prepare the generator - mEventGenerator = Herwig::API::prepareRun(runui); - if (!mEventGenerator) - { - std::cerr << "Error: prepareRun() returned null." << std::endl; - return kFALSE; - } - - mEventGenerator->initialize(); - std::cout << "Herwig generator initialized successfully." << std::endl; - return kTRUE; + using namespace ThePEG; + + if (!std::ifstream(mConfigFile)) + { + LOG(info) << "Run file does not exist: " << mConfigFile; + return kFALSE; } - catch (const std::exception &ex) { - LOG(error) << "Standard exception: " << ex.what(); + SimpleHerwigUI runui(mConfigFile, Herwig::RunMode::RUN, mSeed); + // Prepare the generator + mEventGenerator = Herwig::API::prepareRun(runui); + if (!mEventGenerator) + { + LOG(fatal) << "Error: prepareRun() returned null."; return kFALSE; } + + mEventGenerator->initialize(); + LOG(info) << "Herwig generator initialized successfully."; + return kTRUE; } /// Convert ThePEG event to TParticle format @@ -396,10 +323,11 @@ private: } // namespace eventgen } // namespace o2 -/// Factory function to create HERWIG7 generator from .in file -/// @param inputFile HERWIG input/run file (e.g., "LHC.in"/"LHC.run") -FairGenerator* generateHerwig7(const std::string& inputFile = "LHC.in") +/// HERWIG7 generator from .in/.run file. If seed is -1, a random seed is chosen +/// based on the SimConfig starting seed. +FairGenerator* generateHerwig7(const std::string inputFile = "LHC.in", int seed = -1) { - auto generator = new o2::eventgen::GeneratorHerwig(inputFile); + auto filePath = gSystem->ExpandPathName(inputFile.c_str()); + auto generator = new o2::eventgen::GeneratorHerwig(filePath, seed); return generator; } \ No newline at end of file diff --git a/MC/config/examples/ini/GeneratorHERWIG7.ini b/MC/config/examples/ini/GeneratorHERWIG7.ini index 1d31652fc..af904d1e8 100644 --- a/MC/config/examples/ini/GeneratorHERWIG7.ini +++ b/MC/config/examples/ini/GeneratorHERWIG7.ini @@ -1,3 +1,3 @@ [GeneratorExternal] -fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/examples/epos4herwig/generator_Herwig.C -funcName=generateHerwig7() \ No newline at end of file +fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/examples/herwig/generator_Herwig.C +funcName=generateHerwig7("${O2DPG_MC_CONFIG_ROOT}/MC/config/examples/herwig/LHC.in") diff --git a/MC/config/examples/ini/tests/GeneratorHERWIG7.C b/MC/config/examples/ini/tests/GeneratorHERWIG7.C new file mode 100644 index 000000000..384dda574 --- /dev/null +++ b/MC/config/examples/ini/tests/GeneratorHERWIG7.C @@ -0,0 +1,38 @@ +int External() +{ + std::string path{"o2sim_Kine.root"}; + // Check that file exists, can be opened and has the correct tree + TFile file(path.c_str(), "READ"); + if (file.IsZombie()) + { + std::cerr << "Cannot open ROOT file " << path << "\n"; + return 1; + } + auto tree = (TTree *)file.Get("o2sim"); + if (!tree) + { + std::cerr << "Cannot find tree o2sim in file " << path << "\n"; + return 1; + } + std::vector *tracks{}; + tree->SetBranchAddress("MCTrack", &tracks); + + // Check if all events are filled + auto nEvents = tree->GetEntries(); + for (Long64_t i = 0; i < nEvents; ++i) + { + tree->GetEntry(i); + if (tracks->empty()) + { + std::cerr << "Empty entry found at event " << i << "\n"; + return 1; + } + } + // Check if there are 100 events, as simulated in the o2dpg-test + if (nEvents != 100) + { + std::cerr << "Expected 100 events, got " << nEvents << "\n"; + return 1; + } + return 0; +} \ No newline at end of file