File indexing completed on 2023-09-05 23:29:49
0001 #include <memory>
0002 #include <sstream>
0003 #include <fstream>
0004
0005 #include <HepMC/GenEvent.h>
0006 #include <HepMC/IO_BaseClass.h>
0007
0008 #include <ThePEG/Repository/Repository.h>
0009 #include <ThePEG/EventRecord/Event.h>
0010 #include <ThePEG/Config/ThePEG.h>
0011 #include <ThePEG/LesHouches/LesHouchesReader.h>
0012
0013 #include "FWCore/Framework/interface/LuminosityBlock.h"
0014 #include "FWCore/Framework/interface/EventSetup.h"
0015 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0016
0017 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
0018 #include "SimDataFormats/GeneratorProducts/interface/GenEventInfoProduct.h"
0019 #include "SimDataFormats/GeneratorProducts/interface/GenRunInfoProduct.h"
0020
0021 #include "GeneratorInterface/Core/interface/BaseHadronizer.h"
0022 #include "GeneratorInterface/Core/interface/GeneratorFilter.h"
0023 #include "GeneratorInterface/Core/interface/HadronizerFilter.h"
0024
0025 #include "GeneratorInterface/LHEInterface/interface/LHEProxy.h"
0026
0027 #include "GeneratorInterface/Herwig7Interface/interface/Herwig7Interface.h"
0028
0029 #include <Herwig/API/HerwigAPI.h>
0030 #include "CLHEP/Random/RandomEngine.h"
0031
0032 namespace CLHEP {
0033 class HepRandomEngine;
0034 }
0035
0036 class Herwig7Hadronizer : public Herwig7Interface, public gen::BaseHadronizer {
0037 public:
0038 Herwig7Hadronizer(const edm::ParameterSet& params);
0039 ~Herwig7Hadronizer() override;
0040
0041 bool readSettings(int) { return true; }
0042 bool initializeForInternalPartons();
0043 bool initializeForExternalPartons();
0044 bool declareStableParticles(const std::vector<int>& pdgIds);
0045 bool declareSpecialSettings(const std::vector<std::string>) { return true; }
0046
0047 void statistics();
0048
0049 bool generatePartonsAndHadronize();
0050 bool hadronize();
0051 bool decay();
0052 bool residualDecay();
0053 void finalizeEvent();
0054
0055 const char* classname() const { return "Herwig7Hadronizer"; }
0056 std::unique_ptr<GenLumiInfoHeader> getGenLumiInfoHeader() const override;
0057 void beginLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const&);
0058 void randomizeIndex(edm::LuminosityBlock const& lumi, CLHEP::HepRandomEngine* rengine);
0059
0060 private:
0061 void doSetRandomEngine(CLHEP::HepRandomEngine* v) override { setPEGRandomEngine(v); }
0062
0063 unsigned int eventsToPrint;
0064
0065 ThePEG::EventPtr thepegEvent;
0066 bool haveEvt = false;
0067
0068 std::shared_ptr<lhef::LHEProxy> proxy_;
0069 const std::string handlerDirectory_;
0070 edm::ParameterSet paramSettings;
0071 const std::string runFileName;
0072
0073 unsigned int firstLumiBlock = 0;
0074 unsigned int currentLumiBlock = 0;
0075 };
0076
0077 Herwig7Hadronizer::Herwig7Hadronizer(const edm::ParameterSet& pset)
0078 : Herwig7Interface(pset),
0079 BaseHadronizer(pset),
0080 eventsToPrint(pset.getUntrackedParameter<unsigned int>("eventsToPrint", 0)),
0081 handlerDirectory_(pset.getParameter<std::string>("eventHandlers")),
0082 runFileName(pset.getParameter<std::string>("run")) {
0083 initRepository(pset);
0084 paramSettings = pset;
0085 }
0086
0087 Herwig7Hadronizer::~Herwig7Hadronizer() {}
0088
0089 bool Herwig7Hadronizer::initializeForInternalPartons() {
0090 if (currentLumiBlock == firstLumiBlock) {
0091 std::ifstream runFile(runFileName + ".run");
0092 if (runFile.fail())
0093 {
0094 initRepository(paramSettings);
0095 }
0096 if (!initGenerator()) {
0097 edm::LogInfo("Generator|Herwig7Hadronizer") << "No run step for Herwig chosen. Program will be aborted.";
0098 exit(0);
0099 }
0100 }
0101 return true;
0102 }
0103
0104 bool Herwig7Hadronizer::initializeForExternalPartons() {
0105 if (currentLumiBlock == firstLumiBlock) {
0106 std::ifstream runFile(runFileName + ".run");
0107 if (runFile.fail())
0108 {
0109 initRepository(paramSettings);
0110 }
0111 if (!initGenerator()) {
0112 edm::LogInfo("Generator|Herwig7Hadronizer") << "No run step for Herwig chosen. Program will be aborted.";
0113 exit(0);
0114 }
0115 }
0116 return true;
0117 }
0118
0119 bool Herwig7Hadronizer::declareStableParticles(const std::vector<int>& pdgIds) { return false; }
0120
0121 void Herwig7Hadronizer::statistics() {
0122 if (eg_) {
0123 runInfo().setInternalXSec(
0124 GenRunInfoProduct::XSec(eg_->integratedXSec() / ThePEG::picobarn, eg_->integratedXSecErr() / ThePEG::picobarn));
0125 }
0126 }
0127
0128 bool Herwig7Hadronizer::generatePartonsAndHadronize() {
0129 edm::LogInfo("Generator|Herwig7Hadronizer") << "Start production";
0130
0131 try {
0132 thepegEvent = eg_->shoot();
0133 } catch (std::exception& exc) {
0134 edm::LogWarning("Generator|Herwig7Hadronizer")
0135 << "EGPtr::shoot() thrown an exception, event skipped: " << exc.what();
0136 return false;
0137 }
0138
0139 if (!thepegEvent) {
0140 edm::LogWarning("Generator|Herwig7Hadronizer") << "thepegEvent not initialized";
0141 return false;
0142 }
0143
0144 event() = convert(thepegEvent);
0145 if (!event().get()) {
0146 edm::LogWarning("Generator|Herwig7Hadronizer") << "genEvent not initialized";
0147 return false;
0148 }
0149
0150 return true;
0151 }
0152
0153 bool Herwig7Hadronizer::hadronize() {
0154 if (!haveEvt) {
0155 try {
0156 thepegEvent = eg_->shoot();
0157 haveEvt = true;
0158 } catch (std::exception& exc) {
0159 edm::LogWarning("Generator|Herwig7Hadronizer")
0160 << "EGPtr::shoot() thrown an exception, event skipped: " << exc.what();
0161 return false;
0162 }
0163 }
0164 int evtnum = lheEvent()->evtnum();
0165 if (evtnum == -1) {
0166 edm::LogError("Generator|Herwig7Hadronizer")
0167 << "Event number not set in lhe file, needed for correctly aligning Herwig and LHE events!";
0168 return false;
0169 }
0170 if (thepegEvent->number() < evtnum) {
0171 edm::LogError("Herwig7 interface") << "Herwig does not seem to be generating events in order, did you set "
0172 "/Herwig/EventHandlers/FxFxLHReader:AllowedToReOpen Yes?";
0173 return false;
0174 } else if (thepegEvent->number() == evtnum) {
0175 haveEvt = false;
0176 if (!thepegEvent) {
0177 edm::LogWarning("Generator|Herwig7Hadronizer") << "thepegEvent not initialized";
0178 return false;
0179 }
0180
0181 event() = convert(thepegEvent);
0182 if (!event().get()) {
0183 edm::LogWarning("Generator|Herwig7Hadronizer") << "genEvent not initialized";
0184 return false;
0185 }
0186 return true;
0187 }
0188 edm::LogWarning("Generator|Herwig7Hadronizer") << "Event " << evtnum << " not generated (likely skipped in merging)";
0189 return false;
0190 }
0191
0192 void Herwig7Hadronizer::finalizeEvent() {
0193 eventInfo() = std::make_unique<GenEventInfoProduct>(event().get());
0194 eventInfo()->setBinningValues(std::vector<double>(1, pthat(thepegEvent)));
0195
0196 if (eventsToPrint) {
0197 eventsToPrint--;
0198 event()->print();
0199 }
0200
0201 if (iobc_.get())
0202 iobc_->write_event(event().get());
0203
0204 edm::LogInfo("Generator|Herwig7Hadronizer") << "Event produced";
0205 }
0206
0207 bool Herwig7Hadronizer::decay() { return true; }
0208
0209 bool Herwig7Hadronizer::residualDecay() { return true; }
0210
0211 std::unique_ptr<GenLumiInfoHeader> Herwig7Hadronizer::getGenLumiInfoHeader() const {
0212 auto genLumiInfoHeader = BaseHadronizer::getGenLumiInfoHeader();
0213
0214 if (thepegEvent) {
0215 int weights_number = thepegEvent->optionalWeights().size();
0216
0217 if (weights_number > 1) {
0218 genLumiInfoHeader->weightNames().reserve(weights_number + 1);
0219 genLumiInfoHeader->weightNames().push_back("nominal");
0220 std::map<std::string, double> weights_map = thepegEvent->optionalWeights();
0221 for (std::map<std::string, double>::iterator it = weights_map.begin(); it != weights_map.end(); it++) {
0222 genLumiInfoHeader->weightNames().push_back(it->first);
0223 }
0224 }
0225 }
0226
0227 return genLumiInfoHeader;
0228 }
0229
0230 void Herwig7Hadronizer::randomizeIndex(edm::LuminosityBlock const& lumi, CLHEP::HepRandomEngine* rengine) {
0231 BaseHadronizer::randomizeIndex(lumi, rengine);
0232
0233 if (firstLumiBlock == 0) {
0234 firstLumiBlock = lumi.id().luminosityBlock();
0235 }
0236 currentLumiBlock = lumi.id().luminosityBlock();
0237 }
0238
0239 #include "GeneratorInterface/ExternalDecays/interface/ExternalDecayDriver.h"
0240
0241 typedef edm::GeneratorFilter<Herwig7Hadronizer, gen::ExternalDecayDriver> Herwig7GeneratorFilter;
0242 DEFINE_FWK_MODULE(Herwig7GeneratorFilter);
0243
0244 typedef edm::HadronizerFilter<Herwig7Hadronizer, gen::ExternalDecayDriver> Herwig7HadronizerFilter;
0245 DEFINE_FWK_MODULE(Herwig7HadronizerFilter);