Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 13:06:53

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 
0067   std::shared_ptr<lhef::LHEProxy> proxy_;
0068   const std::string handlerDirectory_;
0069   edm::ParameterSet paramSettings;
0070   const std::string runFileName;
0071 
0072   unsigned int firstLumiBlock = 0;
0073   unsigned int currentLumiBlock = 0;
0074 };
0075 
0076 Herwig7Hadronizer::Herwig7Hadronizer(const edm::ParameterSet& pset)
0077     : Herwig7Interface(pset),
0078       BaseHadronizer(pset),
0079       eventsToPrint(pset.getUntrackedParameter<unsigned int>("eventsToPrint", 0)),
0080       handlerDirectory_(pset.getParameter<std::string>("eventHandlers")),
0081       runFileName(pset.getParameter<std::string>("run")) {
0082   initRepository(pset);
0083   paramSettings = pset;
0084 }
0085 
0086 Herwig7Hadronizer::~Herwig7Hadronizer() {}
0087 
0088 bool Herwig7Hadronizer::initializeForInternalPartons() {
0089   if (currentLumiBlock == firstLumiBlock) {
0090     std::ifstream runFile(runFileName + ".run");
0091     if (runFile.fail())  //required for showering of LHE files
0092     {
0093       initRepository(paramSettings);
0094     }
0095     if (!initGenerator()) {
0096       edm::LogInfo("Generator|Herwig7Hadronizer") << "No run step for Herwig chosen. Program will be aborted.";
0097       exit(0);
0098     }
0099   }
0100   return true;
0101 }
0102 
0103 bool Herwig7Hadronizer::initializeForExternalPartons() {
0104   edm::LogError("Herwig7 interface")
0105       << "Read in of LHE files is not supported in this way. You can read them manually if necessary.";
0106   return false;
0107 }
0108 
0109 bool Herwig7Hadronizer::declareStableParticles(const std::vector<int>& pdgIds) { return false; }
0110 
0111 void Herwig7Hadronizer::statistics() {
0112   if (eg_) {
0113     runInfo().setInternalXSec(
0114         GenRunInfoProduct::XSec(eg_->integratedXSec() / ThePEG::picobarn, eg_->integratedXSecErr() / ThePEG::picobarn));
0115   }
0116 }
0117 
0118 bool Herwig7Hadronizer::generatePartonsAndHadronize() {
0119   edm::LogInfo("Generator|Herwig7Hadronizer") << "Start production";
0120 
0121   try {
0122     thepegEvent = eg_->shoot();
0123   } catch (std::exception& exc) {
0124     edm::LogWarning("Generator|Herwig7Hadronizer")
0125         << "EGPtr::shoot() thrown an exception, event skipped: " << exc.what();
0126     return false;
0127   }
0128 
0129   if (!thepegEvent) {
0130     edm::LogWarning("Generator|Herwig7Hadronizer") << "thepegEvent not initialized";
0131     return false;
0132   }
0133 
0134   event() = convert(thepegEvent);
0135   if (!event().get()) {
0136     edm::LogWarning("Generator|Herwig7Hadronizer") << "genEvent not initialized";
0137     return false;
0138   }
0139 
0140   return true;
0141 }
0142 
0143 bool Herwig7Hadronizer::hadronize() {
0144   edm::LogError("Herwig7 interface")
0145       << "Read in of LHE files is not supported in this way. You can read them manually if necessary.";
0146   return false;
0147 }
0148 
0149 void Herwig7Hadronizer::finalizeEvent() {
0150   eventInfo() = std::make_unique<GenEventInfoProduct>(event().get());
0151   eventInfo()->setBinningValues(std::vector<double>(1, pthat(thepegEvent)));
0152 
0153   if (eventsToPrint) {
0154     eventsToPrint--;
0155     event()->print();
0156   }
0157 
0158   if (iobc_.get())
0159     iobc_->write_event(event().get());
0160 
0161   edm::LogInfo("Generator|Herwig7Hadronizer") << "Event produced";
0162 }
0163 
0164 bool Herwig7Hadronizer::decay() { return true; }
0165 
0166 bool Herwig7Hadronizer::residualDecay() { return true; }
0167 
0168 std::unique_ptr<GenLumiInfoHeader> Herwig7Hadronizer::getGenLumiInfoHeader() const {
0169   auto genLumiInfoHeader = BaseHadronizer::getGenLumiInfoHeader();
0170 
0171   if (thepegEvent) {
0172     int weights_number = thepegEvent->optionalWeights().size();
0173 
0174     if (weights_number > 1) {
0175       genLumiInfoHeader->weightNames().reserve(weights_number + 1);
0176       genLumiInfoHeader->weightNames().push_back("nominal");
0177       std::map<std::string, double> weights_map = thepegEvent->optionalWeights();
0178       for (std::map<std::string, double>::iterator it = weights_map.begin(); it != weights_map.end(); it++) {
0179         genLumiInfoHeader->weightNames().push_back(it->first);
0180       }
0181     }
0182   }
0183 
0184   return genLumiInfoHeader;
0185 }
0186 
0187 void Herwig7Hadronizer::randomizeIndex(edm::LuminosityBlock const& lumi, CLHEP::HepRandomEngine* rengine) {
0188   BaseHadronizer::randomizeIndex(lumi, rengine);
0189 
0190   if (firstLumiBlock == 0) {
0191     firstLumiBlock = lumi.id().luminosityBlock();
0192   }
0193   currentLumiBlock = lumi.id().luminosityBlock();
0194 }
0195 
0196 #include "GeneratorInterface/ExternalDecays/interface/ExternalDecayDriver.h"
0197 
0198 typedef edm::GeneratorFilter<Herwig7Hadronizer, gen::ExternalDecayDriver> Herwig7GeneratorFilter;
0199 DEFINE_FWK_MODULE(Herwig7GeneratorFilter);
0200 
0201 typedef edm::HadronizerFilter<Herwig7Hadronizer, gen::ExternalDecayDriver> Herwig7HadronizerFilter;
0202 DEFINE_FWK_MODULE(Herwig7HadronizerFilter);