Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-07-19 00:58:38

0001 #include <iostream>
0002 #include <memory>
0003 
0004 #include <cstdint>
0005 #include <sstream>
0006 #include <string>
0007 #include <vector>
0008 
0009 #include "HepMC3/GenEvent.h"
0010 #include "HepMC3/GenParticle.h"
0011 #include "HepMC3/Print.h"
0012 
0013 #include "Pythia8/Pythia.h"
0014 #include "Pythia8Plugins/HepMC3.h"
0015 
0016 using namespace Pythia8;
0017 
0018 #include "GeneratorInterface/Pythia8Interface/interface/Py8HMC3InterfaceBase.h"
0019 
0020 #include "GeneratorInterface/Pythia8Interface/plugins/ReweightUserHooks.h"
0021 #include "GeneratorInterface/Pythia8Interface/interface/CustomHook.h"
0022 
0023 // PS matchning prototype
0024 //
0025 #include "GeneratorInterface/Pythia8Interface/plugins/JetMatchingHook.h"
0026 #include "Pythia8Plugins/JetMatching.h"
0027 #include "Pythia8Plugins/aMCatNLOHooks.h"
0028 
0029 // Emission Veto Hooks
0030 //
0031 #include "Pythia8Plugins/PowhegHooks.h"
0032 #include "GeneratorInterface/Pythia8Interface/plugins/EmissionVetoHook1.h"
0033 
0034 // Resonance scale hook
0035 #include "GeneratorInterface/Pythia8Interface/plugins/PowhegResHook.h"
0036 #include "GeneratorInterface/Pythia8Interface/plugins/PowhegHooksBB4L.h"
0037 
0038 //biased tau decayer
0039 #include "GeneratorInterface/Pythia8Interface/interface/BiasedTauDecayer.h"
0040 
0041 //decay filter hook
0042 #include "GeneratorInterface/Pythia8Interface/interface/ResonanceDecayFilterHook.h"
0043 
0044 //decay filter hook
0045 #include "GeneratorInterface/Pythia8Interface/interface/PTFilterHook.h"
0046 
0047 // EvtGen plugin
0048 //
0049 #include "Pythia8Plugins/EvtGen.h"
0050 
0051 #include "FWCore/Concurrency/interface/SharedResourceNames.h"
0052 #include "FWCore/ServiceRegistry/interface/Service.h"
0053 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0054 #include "FWCore/Utilities/interface/RandomNumberGenerator.h"
0055 #include "FWCore/ParameterSet/interface/FileInPath.h"
0056 
0057 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
0058 #include "SimDataFormats/GeneratorProducts/interface/GenRunInfoProduct.h"
0059 
0060 #include "GeneratorInterface/Core/interface/GeneratorFilter.h"
0061 #include "GeneratorInterface/Core/interface/HadronizerFilter.h"
0062 #include "GeneratorInterface/Core/interface/ConcurrentGeneratorFilter.h"
0063 #include "GeneratorInterface/Core/interface/ConcurrentHadronizerFilter.h"
0064 
0065 #include "GeneratorInterface/Pythia8Interface/plugins/LHAupLesHouches.h"
0066 
0067 #include "HepPID/ParticleIDTranslations.hh"
0068 
0069 #include "GeneratorInterface/ExternalDecays/interface/ExternalDecayDriver.h"
0070 #include "GeneratorInterface/ExternalDecays/interface/ConcurrentExternalDecayDriver.h"
0071 
0072 namespace CLHEP {
0073   class HepRandomEngine;
0074 }
0075 
0076 using namespace gen;
0077 
0078 class Pythia8HepMC3Hadronizer : public Py8HMC3InterfaceBase {
0079 public:
0080   Pythia8HepMC3Hadronizer(const edm::ParameterSet &params);
0081   ~Pythia8HepMC3Hadronizer() override = default;
0082 
0083   bool initializeForInternalPartons() override;
0084   bool initializeForExternalPartons();
0085 
0086   bool generatePartonsAndHadronize() override;
0087   bool hadronize();
0088 
0089   virtual bool residualDecay();
0090 
0091   void finalizeEvent() override;
0092 
0093   void statistics() override;
0094 
0095   const char *classname() const override { return "Pythia8HepMC3Hadronizer"; }
0096 
0097   std::unique_ptr<GenLumiInfoHeader> getGenLumiInfoHeader() const override;
0098 
0099 private:
0100   void doSetRandomEngine(CLHEP::HepRandomEngine *v) override { p8SetRandomEngine(v); }
0101   std::vector<std::string> const &doSharedResources() const override { return p8SharedResources; }
0102 
0103   /// Center-of-Mass energy
0104   double comEnergy;
0105 
0106   std::string LHEInputFileName;
0107   std::shared_ptr<LHAupLesHouches> lhaUP;
0108 
0109   enum { PP, PPbar, ElectronPositron };
0110   int fInitialState;  // pp, ppbar, or e-e+
0111 
0112   double fBeam1PZ;
0113   double fBeam2PZ;
0114 
0115   //helper class to allow multiple user hooks simultaneously
0116   std::shared_ptr<UserHooksVector> fUserHooksVector;
0117   bool UserHooksSet;
0118 
0119   // Reweight user hooks
0120   //
0121   std::shared_ptr<UserHooks> fReweightUserHook;
0122   std::shared_ptr<UserHooks> fReweightEmpUserHook;
0123   std::shared_ptr<UserHooks> fReweightRapUserHook;
0124   std::shared_ptr<UserHooks> fReweightPtHatRapUserHook;
0125 
0126   // PS matching prototype
0127   //
0128   std::shared_ptr<JetMatchingHook> fJetMatchingHook;
0129   std::shared_ptr<Pythia8::JetMatchingMadgraph> fJetMatchingPy8InternalHook;
0130   std::shared_ptr<Pythia8::amcnlo_unitarised_interface> fMergingHook;
0131 
0132   // Emission Veto Hooks
0133   //
0134   std::shared_ptr<PowhegHooks> fEmissionVetoHook;
0135   std::shared_ptr<EmissionVetoHook1> fEmissionVetoHook1;
0136 
0137   // Resonance scale hook
0138   std::shared_ptr<PowhegResHook> fPowhegResHook;
0139   std::shared_ptr<PowhegHooksBB4L> fPowhegHooksBB4L;
0140 
0141   // biased tau decayer
0142   std::shared_ptr<BiasedTauDecayer> fBiasedTauDecayer;
0143 
0144   //resonance decay filter hook
0145   std::shared_ptr<ResonanceDecayFilterHook> fResonanceDecayFilterHook;
0146 
0147   //PT filter hook
0148   std::shared_ptr<PTFilterHook> fPTFilterHook;
0149 
0150   //Generic customized hooks vector
0151   std::shared_ptr<UserHooksVector> fCustomHooksVector;
0152 
0153   int EV1_nFinal;
0154   bool EV1_vetoOn;
0155   int EV1_maxVetoCount;
0156   int EV1_pThardMode;
0157   int EV1_pTempMode;
0158   int EV1_emittedMode;
0159   int EV1_pTdefMode;
0160   bool EV1_MPIvetoOn;
0161   int EV1_QEDvetoMode;
0162   int EV1_nFinalMode;
0163 
0164   static const std::vector<std::string> p8SharedResources;
0165 
0166   vector<float> DJR;
0167   int nME;
0168   int nMEFiltered;
0169 
0170   int nISRveto;
0171   int nFSRveto;
0172 };
0173 
0174 const std::vector<std::string> Pythia8HepMC3Hadronizer::p8SharedResources = {edm::SharedResourceNames::kPythia8};
0175 
0176 Pythia8HepMC3Hadronizer::Pythia8HepMC3Hadronizer(const edm::ParameterSet &params)
0177     : Py8HMC3InterfaceBase(params),
0178       comEnergy(params.getParameter<double>("comEnergy")),
0179       LHEInputFileName(params.getUntrackedParameter<std::string>("LHEInputFileName", "")),
0180       fInitialState(PP),
0181       UserHooksSet(false),
0182       nME(-1),
0183       nMEFiltered(-1),
0184       nISRveto(0),
0185       nFSRveto(0) {
0186   ivhepmc = 3;
0187   // J.Y.: the following 3 parameters are hacked "for a reason"
0188   //
0189   if (params.exists("PPbarInitialState")) {
0190     if (fInitialState == PP) {
0191       fInitialState = PPbar;
0192       edm::LogImportant("GeneratorInterface|Pythia8Interface")
0193           << "Pythia8 will be initialized for PROTON-ANTIPROTON INITIAL STATE. "
0194           << "This is a user-request change from the DEFAULT PROTON-PROTON initial state.";
0195     } else {
0196       // probably need to throw on attempt to override ?
0197     }
0198   } else if (params.exists("ElectronPositronInitialState")) {
0199     if (fInitialState == PP) {
0200       fInitialState = ElectronPositron;
0201       edm::LogInfo("GeneratorInterface|Pythia8Interface")
0202           << "Pythia8 will be initialized for ELECTRON-POSITRON INITIAL STATE. "
0203           << "This is a user-request change from the DEFAULT PROTON-PROTON initial state.";
0204     } else {
0205       // probably need to throw on attempt to override ?
0206     }
0207   } else if (params.exists("ElectronProtonInitialState") || params.exists("PositronProtonInitialState")) {
0208     // throw on unknown initial state !
0209     throw edm::Exception(edm::errors::Configuration, "Pythia8Interface")
0210         << " UNKNOWN INITIAL STATE. \n The allowed initial states are: PP, PPbar, ElectronPositron \n";
0211   }
0212 
0213   // avoid filling weights twice (from v8.30x)
0214   toHepMC.set_store_weights(false);
0215 
0216   // Reweight user hook
0217   //
0218   if (params.exists("reweightGen")) {
0219     edm::LogInfo("Pythia8Interface") << "Start setup for reweightGen";
0220     edm::ParameterSet rgParams = params.getParameter<edm::ParameterSet>("reweightGen");
0221     fReweightUserHook.reset(
0222         new PtHatReweightUserHook(rgParams.getParameter<double>("pTRef"), rgParams.getParameter<double>("power")));
0223     edm::LogInfo("Pythia8Interface") << "End setup for reweightGen";
0224   }
0225   if (params.exists("reweightGenEmp")) {
0226     edm::LogInfo("Pythia8Interface") << "Start setup for reweightGenEmp";
0227     edm::ParameterSet rgeParams = params.getParameter<edm::ParameterSet>("reweightGenEmp");
0228 
0229     std::string tuneName = "";
0230     if (rgeParams.exists("tune"))
0231       tuneName = rgeParams.getParameter<std::string>("tune");
0232     fReweightEmpUserHook.reset(new PtHatEmpReweightUserHook(tuneName));
0233     edm::LogInfo("Pythia8Interface") << "End setup for reweightGenEmp";
0234   }
0235   if (params.exists("reweightGenRap")) {
0236     edm::LogInfo("Pythia8Interface") << "Start setup for reweightGenRap";
0237     edm::ParameterSet rgrParams = params.getParameter<edm::ParameterSet>("reweightGenRap");
0238     fReweightRapUserHook.reset(new RapReweightUserHook(rgrParams.getParameter<std::string>("yLabSigmaFunc"),
0239                                                        rgrParams.getParameter<double>("yLabPower"),
0240                                                        rgrParams.getParameter<std::string>("yCMSigmaFunc"),
0241                                                        rgrParams.getParameter<double>("yCMPower"),
0242                                                        rgrParams.getParameter<double>("pTHatMin"),
0243                                                        rgrParams.getParameter<double>("pTHatMax")));
0244     edm::LogInfo("Pythia8Interface") << "End setup for reweightGenRap";
0245   }
0246   if (params.exists("reweightGenPtHatRap")) {
0247     edm::LogInfo("Pythia8Interface") << "Start setup for reweightGenPtHatRap";
0248     edm::ParameterSet rgrParams = params.getParameter<edm::ParameterSet>("reweightGenPtHatRap");
0249     fReweightPtHatRapUserHook.reset(new PtHatRapReweightUserHook(rgrParams.getParameter<std::string>("yLabSigmaFunc"),
0250                                                                  rgrParams.getParameter<double>("yLabPower"),
0251                                                                  rgrParams.getParameter<std::string>("yCMSigmaFunc"),
0252                                                                  rgrParams.getParameter<double>("yCMPower"),
0253                                                                  rgrParams.getParameter<double>("pTHatMin"),
0254                                                                  rgrParams.getParameter<double>("pTHatMax")));
0255     edm::LogInfo("Pythia8Interface") << "End setup for reweightGenPtHatRap";
0256   }
0257 
0258   if (params.exists("useUserHook"))
0259     throw edm::Exception(edm::errors::Configuration, "Pythia8Interface")
0260         << " Obsolete parameter: useUserHook \n Please use the actual one instead \n";
0261 
0262   // PS matching prototype
0263   //
0264   if (params.exists("jetMatching")) {
0265     edm::ParameterSet jmParams = params.getUntrackedParameter<edm::ParameterSet>("jetMatching");
0266     std::string scheme = jmParams.getParameter<std::string>("scheme");
0267     if (scheme == "Madgraph" || scheme == "MadgraphFastJet") {
0268       fJetMatchingHook.reset(new JetMatchingHook(jmParams, &fMasterGen->info));
0269     }
0270   }
0271 
0272   // Pythia8Interface emission veto
0273   //
0274   if (params.exists("emissionVeto1")) {
0275     EV1_nFinal = -1;
0276     if (params.exists("EV1_nFinal"))
0277       EV1_nFinal = params.getParameter<int>("EV1_nFinal");
0278     EV1_vetoOn = true;
0279     if (params.exists("EV1_vetoOn"))
0280       EV1_vetoOn = params.getParameter<bool>("EV1_vetoOn");
0281     EV1_maxVetoCount = 10;
0282     if (params.exists("EV1_maxVetoCount"))
0283       EV1_maxVetoCount = params.getParameter<int>("EV1_maxVetoCount");
0284     EV1_pThardMode = 1;
0285     if (params.exists("EV1_pThardMode"))
0286       EV1_pThardMode = params.getParameter<int>("EV1_pThardMode");
0287     EV1_pTempMode = 0;
0288     if (params.exists("EV1_pTempMode"))
0289       EV1_pTempMode = params.getParameter<int>("EV1_pTempMode");
0290     if (EV1_pTempMode > 2 || EV1_pTempMode < 0)
0291       throw edm::Exception(edm::errors::Configuration, "Pythia8Interface") << " Wrong value for EV1_pTempMode code\n";
0292     EV1_emittedMode = 0;
0293     if (params.exists("EV1_emittedMode"))
0294       EV1_emittedMode = params.getParameter<int>("EV1_emittedMode");
0295     EV1_pTdefMode = 1;
0296     if (params.exists("EV1_pTdefMode"))
0297       EV1_pTdefMode = params.getParameter<int>("EV1_pTdefMode");
0298     EV1_MPIvetoOn = false;
0299     if (params.exists("EV1_MPIvetoOn"))
0300       EV1_MPIvetoOn = params.getParameter<bool>("EV1_MPIvetoOn");
0301     EV1_QEDvetoMode = 0;
0302     if (params.exists("EV1_QEDvetoMode"))
0303       EV1_QEDvetoMode = params.getParameter<int>("EV1_QEDvetoMode");
0304     EV1_nFinalMode = 0;
0305     if (params.exists("EV1_nFinalMode"))
0306       EV1_nFinalMode = params.getParameter<int>("EV1_nFinalMode");
0307     fEmissionVetoHook1.reset(new EmissionVetoHook1(EV1_nFinal,
0308                                                    EV1_vetoOn,
0309                                                    EV1_maxVetoCount,
0310                                                    EV1_pThardMode,
0311                                                    EV1_pTempMode,
0312                                                    EV1_emittedMode,
0313                                                    EV1_pTdefMode,
0314                                                    EV1_MPIvetoOn,
0315                                                    EV1_QEDvetoMode,
0316                                                    EV1_nFinalMode,
0317                                                    0));
0318   }
0319 
0320   if (params.exists("UserCustomization")) {
0321     fCustomHooksVector = std::make_shared<UserHooksVector>();
0322     const std::vector<edm::ParameterSet> userParams =
0323         params.getParameter<std::vector<edm::ParameterSet>>("UserCustomization");
0324     for (const auto &pluginParams : userParams) {
0325       (fCustomHooksVector->hooks)
0326           .push_back(
0327               CustomHookFactory::get()->create(pluginParams.getParameter<std::string>("pluginName"), pluginParams));
0328     }
0329   }
0330 
0331   if (params.exists("VinciaPlugin")) {
0332     throw edm::Exception(edm::errors::Configuration, "Pythia8Interface")
0333         << " Obsolete parameter: VinciaPlugin \n Please use the parameter PartonShowers:model instead \n";
0334   }
0335   if (params.exists("DirePlugin")) {
0336     throw edm::Exception(edm::errors::Configuration, "Pythia8Interface")
0337         << " Obsolete parameter: DirePlugin \n Please use the parameter PartonShowers:model instead \n";
0338   }
0339 }
0340 
0341 bool Pythia8HepMC3Hadronizer::initializeForInternalPartons() {
0342   bool status = false, status1 = false;
0343 
0344   if (lheFile_.empty()) {
0345     if (fInitialState == PP)  // default
0346     {
0347       fMasterGen->settings.mode("Beams:idA", 2212);
0348       fMasterGen->settings.mode("Beams:idB", 2212);
0349     } else if (fInitialState == PPbar) {
0350       fMasterGen->settings.mode("Beams:idA", 2212);
0351       fMasterGen->settings.mode("Beams:idB", -2212);
0352     } else if (fInitialState == ElectronPositron) {
0353       fMasterGen->settings.mode("Beams:idA", 11);
0354       fMasterGen->settings.mode("Beams:idB", -11);
0355     } else {
0356       // throw on unknown initial state !
0357       throw edm::Exception(edm::errors::Configuration, "Pythia8Interface")
0358           << " UNKNOWN INITIAL STATE. \n The allowed initial states are: PP, PPbar, ElectronPositron \n";
0359     }
0360     fMasterGen->settings.parm("Beams:eCM", comEnergy);
0361   } else {
0362     fMasterGen->settings.mode("Beams:frameType", 4);
0363     fMasterGen->settings.word("Beams:LHEF", lheFile_);
0364   }
0365 
0366   if (!fUserHooksVector.get())
0367     fUserHooksVector.reset(new UserHooksVector);
0368   (fUserHooksVector->hooks).clear();
0369 
0370   if (fReweightUserHook.get())
0371     (fUserHooksVector->hooks).push_back(fReweightUserHook);
0372   if (fReweightEmpUserHook.get())
0373     (fUserHooksVector->hooks).push_back(fReweightEmpUserHook);
0374   if (fReweightRapUserHook.get())
0375     (fUserHooksVector->hooks).push_back(fReweightRapUserHook);
0376   if (fReweightPtHatRapUserHook.get())
0377     (fUserHooksVector->hooks).push_back(fReweightPtHatRapUserHook);
0378   if (fJetMatchingHook.get())
0379     (fUserHooksVector->hooks).push_back(fJetMatchingHook);
0380   if (fEmissionVetoHook1.get()) {
0381     edm::LogInfo("Pythia8Interface") << "Turning on Emission Veto Hook 1 from CMSSW Pythia8Interface";
0382     (fUserHooksVector->hooks).push_back(fEmissionVetoHook1);
0383   }
0384 
0385   if (fMasterGen->settings.mode("POWHEG:veto") > 0 || fMasterGen->settings.mode("POWHEG:MPIveto") > 0) {
0386     if (fJetMatchingHook.get() || fEmissionVetoHook1.get())
0387       throw edm::Exception(edm::errors::Configuration, "Pythia8Interface")
0388           << " Attempt to turn on PowhegHooks by pythia8 settings but there are incompatible hooks on \n Incompatible "
0389              "are : jetMatching, emissionVeto1 \n";
0390 
0391     if (!fEmissionVetoHook.get())
0392       fEmissionVetoHook.reset(new PowhegHooks());
0393 
0394     edm::LogInfo("Pythia8Interface") << "Turning on Emission Veto Hook from pythia8 code";
0395     (fUserHooksVector->hooks).push_back(fEmissionVetoHook);
0396   }
0397 
0398   bool PowhegRes = fMasterGen->settings.flag("POWHEGres:calcScales");
0399   if (PowhegRes) {
0400     edm::LogInfo("Pythia8Interface") << "Turning on resonance scale setting from CMSSW Pythia8Interface";
0401     if (!fPowhegResHook.get())
0402       fPowhegResHook.reset(new PowhegResHook());
0403     (fUserHooksVector->hooks).push_back(fPowhegResHook);
0404   }
0405 
0406   bool PowhegBB4L = fMasterGen->settings.flag("POWHEG:bb4l");
0407   if (PowhegBB4L) {
0408     edm::LogInfo("Pythia8Interface") << "Turning on BB4l hook from CMSSW Pythia8Interface";
0409     if (!fPowhegHooksBB4L.get())
0410       fPowhegHooksBB4L.reset(new PowhegHooksBB4L());
0411     (fUserHooksVector->hooks).push_back(fPowhegHooksBB4L);
0412   }
0413 
0414   //adapted from main89.cc in pythia8 examples
0415   bool internalMatching = fMasterGen->settings.flag("JetMatching:merge");
0416   bool internalMerging = !(fMasterGen->settings.word("Merging:Process") == "void");
0417 
0418   if (internalMatching && internalMerging) {
0419     throw edm::Exception(edm::errors::Configuration, "Pythia8Interface")
0420         << " Only one jet matching/merging scheme can be used at a time. \n";
0421   }
0422 
0423   if (internalMatching) {
0424     if (!fJetMatchingPy8InternalHook.get())
0425       fJetMatchingPy8InternalHook.reset(new Pythia8::JetMatchingMadgraph);
0426     (fUserHooksVector->hooks).push_back(fJetMatchingPy8InternalHook);
0427   }
0428 
0429   if (internalMerging) {
0430     int scheme = (fMasterGen->settings.flag("Merging:doUMEPSTree") || fMasterGen->settings.flag("Merging:doUMEPSSubt"))
0431                      ? 1
0432                      : ((fMasterGen->settings.flag("Merging:doUNLOPSTree") ||
0433                          fMasterGen->settings.flag("Merging:doUNLOPSSubt") ||
0434                          fMasterGen->settings.flag("Merging:doUNLOPSLoop") ||
0435                          fMasterGen->settings.flag("Merging:doUNLOPSSubtNLO"))
0436                             ? 2
0437                             : 0);
0438     if (!fMergingHook.get())
0439       fMergingHook.reset(new Pythia8::amcnlo_unitarised_interface(scheme));
0440     (fUserHooksVector->hooks).push_back(fMergingHook);
0441   }
0442 
0443   bool biasedTauDecayer = fMasterGen->settings.flag("BiasedTauDecayer:filter");
0444   if (biasedTauDecayer) {
0445     if (!fBiasedTauDecayer.get()) {
0446       Pythia8::Info localInfo = fMasterGen->info;
0447       fBiasedTauDecayer.reset(new BiasedTauDecayer(&localInfo, &(fMasterGen->settings)));
0448     }
0449     std::vector<int> handledParticles;
0450     handledParticles.push_back(15);
0451     fMasterGen->setDecayPtr(fBiasedTauDecayer, handledParticles);
0452   }
0453 
0454   bool resonanceDecayFilter = fMasterGen->settings.flag("ResonanceDecayFilter:filter");
0455   if (resonanceDecayFilter) {
0456     fResonanceDecayFilterHook.reset(new ResonanceDecayFilterHook);
0457     (fUserHooksVector->hooks).push_back(fResonanceDecayFilterHook);
0458   }
0459 
0460   bool PTFilter = fMasterGen->settings.flag("PTFilter:filter");
0461   if (PTFilter) {
0462     fPTFilterHook.reset(new PTFilterHook);
0463     (fUserHooksVector->hooks).push_back(fPTFilterHook);
0464   }
0465 
0466   if (!(fUserHooksVector->hooks).empty() && !UserHooksSet) {
0467     for (auto &fUserHook : fUserHooksVector->hooks) {
0468       fMasterGen->addUserHooksPtr(fUserHook);
0469     }
0470     UserHooksSet = true;
0471   }
0472 
0473   if (fCustomHooksVector.get()) {
0474     edm::LogInfo("Pythia8Interface") << "Adding customized user hooks";
0475     for (const auto &fUserHook : fCustomHooksVector->hooks) {
0476       fMasterGen->addUserHooksPtr(fUserHook);
0477     }
0478   }
0479 
0480   edm::LogInfo("Pythia8Interface") << "Initializing MasterGen";
0481   status = fMasterGen->init();
0482 
0483   //clean up temp file
0484   if (!slhafile_.empty()) {
0485     std::remove(slhafile_.c_str());
0486   }
0487 
0488   if (pythiaPylistVerbosity > 10) {
0489     if (pythiaPylistVerbosity == 11 || pythiaPylistVerbosity == 13)
0490       fMasterGen->settings.listAll();
0491     if (pythiaPylistVerbosity == 12 || pythiaPylistVerbosity == 13)
0492       fMasterGen->particleData.listAll();
0493   }
0494 
0495   // init decayer
0496   fDecayer->settings.flag("ProcessLevel:all", false);  // trick
0497   fDecayer->settings.flag("ProcessLevel:resonanceDecays", true);
0498   edm::LogInfo("Pythia8Interface") << "Initializing Decayer";
0499   status1 = fDecayer->init();
0500 
0501   if (useEvtGen) {
0502     edm::LogInfo("Pythia8Hadronizer") << "Creating and initializing pythia8 EvtGen plugin";
0503     if (!evtgenDecays.get()) {
0504       evtgenDecays.reset(new EvtGenDecays(fMasterGen.get(), evtgenDecFile, evtgenPdlFile));
0505       for (unsigned int i = 0; i < evtgenUserFiles.size(); i++)
0506         evtgenDecays->readDecayFile(evtgenUserFiles.at(i));
0507     }
0508   }
0509 
0510   return (status && status1);
0511 }
0512 
0513 bool Pythia8HepMC3Hadronizer::initializeForExternalPartons() {
0514   edm::LogInfo("Pythia8Interface") << "Initializing for external partons";
0515 
0516   bool status = false, status1 = false;
0517 
0518   if (!fUserHooksVector.get())
0519     fUserHooksVector.reset(new UserHooksVector);
0520   (fUserHooksVector->hooks).clear();
0521 
0522   if (fReweightUserHook.get())
0523     (fUserHooksVector->hooks).push_back(fReweightUserHook);
0524   if (fReweightEmpUserHook.get())
0525     (fUserHooksVector->hooks).push_back(fReweightEmpUserHook);
0526   if (fReweightRapUserHook.get())
0527     (fUserHooksVector->hooks).push_back(fReweightRapUserHook);
0528   if (fReweightPtHatRapUserHook.get())
0529     (fUserHooksVector->hooks).push_back(fReweightPtHatRapUserHook);
0530   if (fJetMatchingHook.get())
0531     (fUserHooksVector->hooks).push_back(fJetMatchingHook);
0532   if (fEmissionVetoHook1.get()) {
0533     edm::LogInfo("Pythia8Interface") << "Turning on Emission Veto Hook 1 from CMSSW Pythia8Interface";
0534     (fUserHooksVector->hooks).push_back(fEmissionVetoHook1);
0535   }
0536 
0537   if (fCustomHooksVector.get()) {
0538     edm::LogInfo("Pythia8Interface") << "Adding customized user hook";
0539     for (const auto &fUserHook : fCustomHooksVector->hooks) {
0540       (fUserHooksVector->hooks).push_back(fUserHook);
0541     }
0542   }
0543 
0544   if (fMasterGen->settings.mode("POWHEG:veto") > 0 || fMasterGen->settings.mode("POWHEG:MPIveto") > 0) {
0545     if (fJetMatchingHook.get() || fEmissionVetoHook1.get())
0546       throw edm::Exception(edm::errors::Configuration, "Pythia8Interface")
0547           << " Attempt to turn on PowhegHooks by pythia8 settings but there are incompatible hooks on \n Incompatible "
0548              "are : jetMatching, emissionVeto1 \n";
0549 
0550     if (!fEmissionVetoHook.get())
0551       fEmissionVetoHook.reset(new PowhegHooks());
0552 
0553     edm::LogInfo("Pythia8Interface") << "Turning on Emission Veto Hook from pythia8 code";
0554     (fUserHooksVector->hooks).push_back(fEmissionVetoHook);
0555   }
0556 
0557   bool PowhegRes = fMasterGen->settings.flag("POWHEGres:calcScales");
0558   if (PowhegRes) {
0559     edm::LogInfo("Pythia8Interface") << "Turning on resonance scale setting from CMSSW Pythia8Interface";
0560     if (!fPowhegResHook.get())
0561       fPowhegResHook.reset(new PowhegResHook());
0562     (fUserHooksVector->hooks).push_back(fPowhegResHook);
0563   }
0564 
0565   bool PowhegBB4L = fMasterGen->settings.flag("POWHEG:bb4l");
0566   if (PowhegBB4L) {
0567     edm::LogInfo("Pythia8Interface") << "Turning on BB4l hook from CMSSW Pythia8Interface";
0568     if (!fPowhegHooksBB4L.get())
0569       fPowhegHooksBB4L.reset(new PowhegHooksBB4L());
0570     (fUserHooksVector->hooks).push_back(fPowhegHooksBB4L);
0571   }
0572 
0573   //adapted from main89.cc in pythia8 examples
0574   bool internalMatching = fMasterGen->settings.flag("JetMatching:merge");
0575   bool internalMerging = !(fMasterGen->settings.word("Merging:Process") == "void");
0576 
0577   if (internalMatching && internalMerging) {
0578     throw edm::Exception(edm::errors::Configuration, "Pythia8Interface")
0579         << " Only one jet matching/merging scheme can be used at a time. \n";
0580   }
0581 
0582   if (internalMatching) {
0583     if (!fJetMatchingPy8InternalHook.get())
0584       fJetMatchingPy8InternalHook.reset(new Pythia8::JetMatchingMadgraph);
0585     (fUserHooksVector->hooks).push_back(fJetMatchingPy8InternalHook);
0586   }
0587 
0588   if (internalMerging) {
0589     int scheme = (fMasterGen->settings.flag("Merging:doUMEPSTree") || fMasterGen->settings.flag("Merging:doUMEPSSubt"))
0590                      ? 1
0591                      : ((fMasterGen->settings.flag("Merging:doUNLOPSTree") ||
0592                          fMasterGen->settings.flag("Merging:doUNLOPSSubt") ||
0593                          fMasterGen->settings.flag("Merging:doUNLOPSLoop") ||
0594                          fMasterGen->settings.flag("Merging:doUNLOPSSubtNLO"))
0595                             ? 2
0596                             : 0);
0597     if (!fMergingHook.get())
0598       fMergingHook.reset(new Pythia8::amcnlo_unitarised_interface(scheme));
0599     (fUserHooksVector->hooks).push_back(fMergingHook);
0600   }
0601 
0602   bool biasedTauDecayer = fMasterGen->settings.flag("BiasedTauDecayer:filter");
0603   if (biasedTauDecayer) {
0604     if (!fBiasedTauDecayer.get()) {
0605       Pythia8::Info localInfo = fMasterGen->info;
0606       fBiasedTauDecayer.reset(new BiasedTauDecayer(&localInfo, &(fMasterGen->settings)));
0607     }
0608     std::vector<int> handledParticles;
0609     handledParticles.push_back(15);
0610     fMasterGen->setDecayPtr(fBiasedTauDecayer, handledParticles);
0611   }
0612 
0613   bool resonanceDecayFilter = fMasterGen->settings.flag("ResonanceDecayFilter:filter");
0614   if (resonanceDecayFilter) {
0615     fResonanceDecayFilterHook.reset(new ResonanceDecayFilterHook);
0616     (fUserHooksVector->hooks).push_back(fResonanceDecayFilterHook);
0617   }
0618 
0619   bool PTFilter = fMasterGen->settings.flag("PTFilter:filter");
0620   if (PTFilter) {
0621     fPTFilterHook.reset(new PTFilterHook);
0622     (fUserHooksVector->hooks).push_back(fPTFilterHook);
0623   }
0624 
0625   if (!(fUserHooksVector->hooks).empty() && !UserHooksSet) {
0626     for (auto &fUserHook : fUserHooksVector->hooks) {
0627       fMasterGen->addUserHooksPtr(fUserHook);
0628     }
0629     UserHooksSet = true;
0630   }
0631 
0632   if (!LHEInputFileName.empty()) {
0633     edm::LogInfo("Pythia8Interface") << "Initialize direct pythia8 reading from LHE file " << LHEInputFileName;
0634     edm::LogInfo("Pythia8Interface") << "Some LHE information can be not stored";
0635     fMasterGen->settings.mode("Beams:frameType", 4);
0636     fMasterGen->settings.word("Beams:LHEF", LHEInputFileName);
0637     status = fMasterGen->init();
0638 
0639   } else {
0640     lhaUP.reset(new LHAupLesHouches());
0641     lhaUP->setScalesFromLHEF(fMasterGen->settings.flag("Beams:setProductionScalesFromLHEF"));
0642     lhaUP->loadRunInfo(lheRunInfo());
0643 
0644     if (fJetMatchingHook.get()) {
0645       fJetMatchingHook->init(lheRunInfo());
0646     }
0647 
0648     fMasterGen->settings.mode("Beams:frameType", 5);
0649     fMasterGen->setLHAupPtr(lhaUP);
0650     edm::LogInfo("Pythia8Interface") << "Initializing MasterGen";
0651     status = fMasterGen->init();
0652   }
0653 
0654   //clean up temp file
0655   if (!slhafile_.empty()) {
0656     std::remove(slhafile_.c_str());
0657   }
0658 
0659   if (pythiaPylistVerbosity > 10) {
0660     if (pythiaPylistVerbosity == 11 || pythiaPylistVerbosity == 13)
0661       fMasterGen->settings.listAll();
0662     if (pythiaPylistVerbosity == 12 || pythiaPylistVerbosity == 13)
0663       fMasterGen->particleData.listAll();
0664   }
0665 
0666   // init decayer
0667   fDecayer->settings.flag("ProcessLevel:all", false);  // trick
0668   fDecayer->settings.flag("ProcessLevel:resonanceDecays", true);
0669   edm::LogInfo("Pythia8Interface") << "Initializing Decayer";
0670   status1 = fDecayer->init();
0671 
0672   if (useEvtGen) {
0673     edm::LogInfo("Pythia8Hadronizer") << "Creating and initializing pythia8 EvtGen plugin";
0674     if (!evtgenDecays.get()) {
0675       evtgenDecays.reset(new EvtGenDecays(fMasterGen.get(), evtgenDecFile, evtgenPdlFile));
0676       for (unsigned int i = 0; i < evtgenUserFiles.size(); i++)
0677         evtgenDecays->readDecayFile(evtgenUserFiles.at(i));
0678     }
0679   }
0680 
0681   return (status && status1);
0682 }
0683 
0684 void Pythia8HepMC3Hadronizer::statistics() {
0685   fMasterGen->stat();
0686 
0687   if (fEmissionVetoHook.get()) {
0688     edm::LogPrint("Pythia8Interface") << "\n"
0689                                       << "Number of ISR vetoed = " << nISRveto;
0690     edm::LogPrint("Pythia8Interface") << "Number of FSR vetoed = " << nFSRveto;
0691   }
0692 
0693   double xsec = fMasterGen->info.sigmaGen();  // cross section in mb
0694   xsec *= 1.0e9;                              // translate to pb (CMS/Gen "convention" as of May 2009)
0695   double err = fMasterGen->info.sigmaErr();   // cross section err in mb
0696   err *= 1.0e9;                               // translate to pb (CMS/Gen "convention" as of May 2009)
0697   runInfo().setInternalXSec(GenRunInfoProduct::XSec(xsec, err));
0698 }
0699 
0700 bool Pythia8HepMC3Hadronizer::generatePartonsAndHadronize() {
0701   DJR.resize(0);
0702   nME = -1;
0703   nMEFiltered = -1;
0704 
0705   if (fJetMatchingHook.get()) {
0706     fJetMatchingHook->resetMatchingStatus();
0707     fJetMatchingHook->beforeHadronization(lheEvent());
0708   }
0709 
0710   if (!fMasterGen->next())
0711     return false;
0712 
0713   double mergeweight = fMasterGen.get()->info.mergingWeightNLO();
0714   if (fMergingHook.get()) {
0715     mergeweight *= fMergingHook->getNormFactor();
0716   }
0717 
0718   //protect against 0-weight from ckkw or similar
0719   if (std::abs(mergeweight) == 0.) {
0720     event().reset();
0721     return false;
0722   }
0723 
0724   if (fJetMatchingPy8InternalHook.get()) {
0725     const std::vector<double> djrmatch = fJetMatchingPy8InternalHook->getDJR();
0726     //cap size of djr vector to save storage space (keep only up to first 6 elements)
0727     unsigned int ndjr = std::min(djrmatch.size(), std::vector<double>::size_type(6));
0728     for (unsigned int idjr = 0; idjr < ndjr; ++idjr) {
0729       DJR.push_back(djrmatch[idjr]);
0730     }
0731 
0732     nME = fJetMatchingPy8InternalHook->nMEpartons().first;
0733     nMEFiltered = fJetMatchingPy8InternalHook->nMEpartons().second;
0734   }
0735 
0736   if (evtgenDecays.get())
0737     evtgenDecays->decay();
0738 
0739   //event() = std::make_unique<HepMC::GenEvent>();
0740   event3() = std::make_unique<HepMC3::GenEvent>();
0741   bool py8hepmc = toHepMC.fill_next_event(*(fMasterGen.get()), event3().get());
0742 
0743   if (!py8hepmc) {
0744     return false;
0745   }
0746 
0747   // 0th weight is not filled anymore since v8.30x, pushback manually
0748   event3()->weights().push_back(fMasterGen->info.weight());
0749 
0750   //add ckkw/umeps/unlops merging weight
0751   if (mergeweight != 1.) {
0752     event3()->weights()[0] *= mergeweight;
0753   }
0754 
0755   if (fEmissionVetoHook.get()) {
0756     nISRveto += fEmissionVetoHook->getNISRveto();
0757     nFSRveto += fEmissionVetoHook->getNFSRveto();
0758   }
0759 
0760   //fill additional weights for systematic uncertainties
0761   if (fMasterGen->info.getWeightsDetailedSize() > 0) {
0762     for (const string &key : fMasterGen->info.initrwgt->weightsKeys) {
0763       double wgt = (*fMasterGen->info.weights_detailed)[key];
0764       event3()->weights().push_back(wgt);
0765     }
0766   } else if (fMasterGen->info.getWeightsCompressedSize() > 0) {
0767     for (unsigned int i = 0; i < fMasterGen->info.getWeightsCompressedSize(); i++) {
0768       double wgt = fMasterGen->info.getWeightsCompressedValue(i);
0769       event3()->weights().push_back(wgt);
0770     }
0771   }
0772 
0773   // fill shower weights
0774   // http://home.thep.lu.se/~torbjorn/pythia82html/Variations.html
0775   if (fMasterGen->info.nWeights() > 1) {
0776     for (int i = 0; i < fMasterGen->info.nWeights(); ++i) {
0777       double wgt = fMasterGen->info.weight(i);
0778       event3()->weights().push_back(wgt);
0779     }
0780   }
0781 
0782 #if 0
0783   // VINCIA shower weights
0784   // http://vincia.hepforge.org/current/share/Vincia/htmldoc/VinciaUncertainties.html
0785   if (fvincia.get()) {
0786     event()->weights()[0] *= fvincia->weight(0);
0787     for (int iVar = 1; iVar < fvincia->nWeights(); iVar++) {
0788       event()->weights().push_back(fvincia->weight(iVar));
0789     }
0790   }
0791 
0792   // Retrieve Dire shower weights
0793   if (fDire.get()) {
0794     fDire->weightsPtr->calcWeight(0.);
0795     fDire->weightsPtr->reset();
0796 
0797     //Make sure the base weight comes first
0798     event()->weights()[0] *= fDire->weightsPtr->getShowerWeight("base");
0799 
0800     unordered_map<string, double>::iterator it;
0801     for (it = fDire->weightsPtr->getShowerWeights()->begin(); it != fDire->weightsPtr->getShowerWeights()->end();
0802          it++) {
0803       if (it->first == "base")
0804         continue;
0805       event()->weights().push_back(it->second);
0806     }
0807   }
0808 #endif
0809 
0810   return true;
0811 }
0812 
0813 bool Pythia8HepMC3Hadronizer::hadronize() {
0814   DJR.resize(0);
0815   nME = -1;
0816   nMEFiltered = -1;
0817   if (LHEInputFileName.empty())
0818     lhaUP->loadEvent(lheEvent());
0819 
0820   if (fJetMatchingHook.get()) {
0821     fJetMatchingHook->resetMatchingStatus();
0822     fJetMatchingHook->beforeHadronization(lheEvent());
0823   }
0824 
0825   bool py8next = fMasterGen->next();
0826 
0827   double mergeweight = fMasterGen.get()->info.mergingWeightNLO();
0828   if (fMergingHook.get()) {
0829     mergeweight *= fMergingHook->getNormFactor();
0830   }
0831 
0832   //protect against 0-weight from ckkw or similar
0833   if (!py8next || std::abs(mergeweight) == 0.) {
0834     lheEvent()->count(lhef::LHERunInfo::kSelected, 1.0, mergeweight);
0835     event3().reset();
0836     return false;
0837   }
0838 
0839   if (fJetMatchingPy8InternalHook.get()) {
0840     const std::vector<double> djrmatch = fJetMatchingPy8InternalHook->getDJR();
0841     //cap size of djr vector to save storage space (keep only up to first 6 elements)
0842     unsigned int ndjr = std::min(djrmatch.size(), std::vector<double>::size_type(6));
0843     for (unsigned int idjr = 0; idjr < ndjr; ++idjr) {
0844       DJR.push_back(djrmatch[idjr]);
0845     }
0846 
0847     nME = fJetMatchingPy8InternalHook->nMEpartons().first;
0848     nMEFiltered = fJetMatchingPy8InternalHook->nMEpartons().second;
0849   }
0850 
0851   // update LHE matching statistics
0852   //
0853   lheEvent()->count(lhef::LHERunInfo::kAccepted, 1.0, mergeweight);
0854 
0855   if (evtgenDecays.get())
0856     evtgenDecays->decay();
0857 
0858   event() = std::make_unique<HepMC::GenEvent>();
0859   event3() = std::make_unique<HepMC3::GenEvent>();
0860   bool py8hepmc = toHepMC.fill_next_event(*(fMasterGen.get()), event3().get());
0861 
0862   if (!py8hepmc) {
0863     return false;
0864   }
0865 
0866   // 0th weight is not filled anymore since v8.30x, pushback manually
0867   event3()->weights().push_back(fMasterGen->info.weight());
0868 
0869   //add ckkw/umeps/unlops merging weight
0870   if (mergeweight != 1.) {
0871     event3()->weights()[0] *= mergeweight;
0872   }
0873 
0874   if (fEmissionVetoHook.get()) {
0875     nISRveto += fEmissionVetoHook->getNISRveto();
0876     nFSRveto += fEmissionVetoHook->getNFSRveto();
0877   }
0878 
0879   // fill shower weights
0880   // http://home.thep.lu.se/~torbjorn/pythia82html/Variations.html
0881   if (fMasterGen->info.nWeights() > 1) {
0882     for (int i = 0; i < fMasterGen->info.nWeights(); ++i) {
0883       double wgt = fMasterGen->info.weight(i);
0884       event3()->weights().push_back(wgt);
0885     }
0886   }
0887 
0888   return true;
0889 }
0890 
0891 bool Pythia8HepMC3Hadronizer::residualDecay() {
0892   Event *pythiaEvent = &(fMasterGen->event);
0893 
0894   int NPartsBeforeDecays = pythiaEvent->size() - 1;  // do NOT count the very 1st "system" particle
0895                                                      // in Pythia8::Event record; it does NOT even
0896                                                      // get translated by the HepMCInterface to the
0897                                                      // HepMC::GenEvent record !!!
0898   //int NPartsAfterDecays = event().get()->particles_size();
0899   int NPartsAfterDecays = 0;
0900   for (auto p : (event3().get())->particles()) {
0901     NPartsAfterDecays++;
0902   }
0903 
0904   if (NPartsAfterDecays == NPartsBeforeDecays)
0905     return true;
0906 
0907   bool result = true;
0908 
0909 // this below is to be rewritten in future development, then it will be uncommented
0910 #if 0
0911   for (int ipart = NPartsAfterDecays; ipart > NPartsBeforeDecays; ipart--) {
0912     HepMC::GenParticle *part = event().get()->barcode_to_particle(ipart);
0913 
0914     if (part->status() == 1 && (fDecayer->particleData).canDecay(part->pdg_id())) {
0915       fDecayer->event.reset();
0916       Particle py8part(part->pdg_id(),
0917                        93,
0918                        0,
0919                        0,
0920                        0,
0921                        0,
0922                        0,
0923                        0,
0924                        part->momentum().x(),
0925                        part->momentum().y(),
0926                        part->momentum().z(),
0927                        part->momentum().t(),
0928                        part->generated_mass());
0929       HepMC::GenVertex *ProdVtx = part->production_vertex();
0930       py8part.vProd(ProdVtx->position().x(), ProdVtx->position().y(), ProdVtx->position().z(), ProdVtx->position().t());
0931       py8part.tau((fDecayer->particleData).tau0(part->pdg_id()));
0932       fDecayer->event.append(py8part);
0933       int nentries = fDecayer->event.size();
0934       if (!fDecayer->event[nentries - 1].mayDecay())
0935         continue;
0936       fDecayer->next();
0937       int nentries1 = fDecayer->event.size();
0938       if (nentries1 <= nentries)
0939         continue;  //same number of particles, no decays...
0940 
0941       part->set_status(2);
0942 
0943       result = toHepMC.fill_next_event(*(fDecayer.get()), event().get(), -1, true, part);
0944     }
0945   }
0946 #endif
0947 
0948   return result;
0949 }
0950 
0951 void Pythia8HepMC3Hadronizer::finalizeEvent() {
0952   bool lhe = lheEvent() != nullptr;
0953 
0954   // protection against empty weight container
0955   if ((event3()->weights()).empty())
0956     (event3()->weights()).push_back(1.);
0957 
0958   // now create the GenEventInfo product from the GenEvent and fill
0959   // the missing pieces
0960   eventInfo3() = std::make_unique<GenEventInfoProduct3>(event3().get());
0961 
0962   // in pythia pthat is used to subdivide samples into different bins
0963   // in LHE mode the binning is done by the external ME generator
0964   // which is likely not pthat, so only filling it for Py6 internal mode
0965   if (!lhe) {
0966     eventInfo3()->setBinningValues(std::vector<double>(1, fMasterGen->info.pTHat()));
0967   }
0968 
0969   eventInfo3()->setDJR(DJR);
0970   eventInfo3()->setNMEPartons(nME);
0971   eventInfo3()->setNMEPartonsFiltered(nMEFiltered);
0972 
0973   //******** Verbosity ********
0974 
0975   if (maxEventsToPrint > 0 && (pythiaPylistVerbosity || pythiaHepMCVerbosity || pythiaHepMCVerbosityParticles)) {
0976     maxEventsToPrint--;
0977     if (pythiaPylistVerbosity) {
0978       fMasterGen->info.list();
0979       fMasterGen->event.list();
0980     }
0981 
0982     if (pythiaHepMCVerbosity) {
0983       std::cout << "Event process = " << fMasterGen->info.code() << "\n"
0984                 << "----------------------" << std::endl;
0985       //event3()->print();
0986       HepMC3::Print::listing(*(event3().get()));
0987     }
0988     if (pythiaHepMCVerbosityParticles) {
0989       std::cout << "Event process = " << fMasterGen->info.code() << "\n"
0990                 << "----------------------" << std::endl;
0991       //ascii_io->write_event(event().get());
0992       for (const auto &p : (event3().get())->particles()) {
0993         HepMC3::Print::line(p, true);
0994       }
0995     }
0996   }
0997 }
0998 
0999 std::unique_ptr<GenLumiInfoHeader> Pythia8HepMC3Hadronizer::getGenLumiInfoHeader() const {
1000   auto genLumiInfoHeader = BaseHadronizer::getGenLumiInfoHeader();
1001 
1002   //fill lhe headers
1003   //*FIXME* initrwgt header is corrupt due to pythia bug
1004   for (const std::string &key : fMasterGen->info.headerKeys()) {
1005     genLumiInfoHeader->lheHeaders().emplace_back(key, fMasterGen->info.header(key));
1006   }
1007 
1008   //check, if it is not only nominal weight
1009   int weights_number = fMasterGen->info.nWeights();
1010   if (fMasterGen->info.initrwgt)
1011     weights_number += fMasterGen->info.initrwgt->weightsKeys.size();
1012   if (weights_number > 1) {
1013     genLumiInfoHeader->weightNames().reserve(weights_number + 1);
1014     genLumiInfoHeader->weightNames().push_back("nominal");
1015   }
1016 
1017   //fill weight names
1018   if (fMasterGen->info.initrwgt) {
1019     for (const std::string &key : fMasterGen->info.initrwgt->weightsKeys) {
1020       std::string weightgroupname;
1021       for (const auto &wgtgrp : fMasterGen->info.initrwgt->weightgroups) {
1022         const auto &wgtgrpwgt = wgtgrp.second.weights.find(key);
1023         if (wgtgrpwgt != wgtgrp.second.weights.end()) {
1024           weightgroupname = wgtgrp.first;
1025         }
1026       }
1027 
1028       std::ostringstream weightname;
1029       weightname << "LHE, id = " << key << ", ";
1030       if (!weightgroupname.empty()) {
1031         weightname << "group = " << weightgroupname << ", ";
1032       }
1033       weightname << fMasterGen->info.initrwgt->weights[key].contents;
1034       genLumiInfoHeader->weightNames().push_back(weightname.str());
1035     }
1036   }
1037 
1038   //fill shower labels
1039   // http://home.thep.lu.se/~torbjorn/pythia82html/Variations.html
1040   // http://home.thep.lu.se/~torbjorn/doxygen/classPythia8_1_1Info.html
1041   if (fMasterGen->info.nWeights() > 1) {
1042     for (int i = 0; i < fMasterGen->info.nWeights(); ++i) {
1043       genLumiInfoHeader->weightNames().push_back(fMasterGen->info.weightLabel(i));
1044     }
1045   }
1046 
1047 #if 0
1048   // VINCIA shower weights
1049   // http://vincia.hepforge.org/current/share/Vincia/htmldoc/VinciaUncertainties.html
1050   if (fvincia.get()) {
1051     for (int iVar = 0; iVar < fvincia->nWeights(); iVar++) {
1052       genLumiInfoHeader->weightNames().push_back(fvincia->weightLabel(iVar));
1053     }
1054   }
1055 
1056   if (fDire.get()) {
1057     //Make sure the base weight comes first
1058     genLumiInfoHeader->weightNames().push_back("base");
1059 
1060     unordered_map<string, double>::iterator it;
1061     for (it = fDire->weightsPtr->getShowerWeights()->begin(); it != fDire->weightsPtr->getShowerWeights()->end();
1062          it++) {
1063       if (it->first == "base")
1064         continue;
1065       genLumiInfoHeader->weightNames().push_back(it->first);
1066     }
1067   }
1068 #endif
1069 
1070   return genLumiInfoHeader;
1071 }
1072 
1073 typedef edm::GeneratorFilter<Pythia8HepMC3Hadronizer, ExternalDecayDriver> Pythia8HepMC3GeneratorFilter;
1074 DEFINE_FWK_MODULE(Pythia8HepMC3GeneratorFilter);
1075 
1076 typedef edm::HadronizerFilter<Pythia8HepMC3Hadronizer, ExternalDecayDriver> Pythia8HepMC3HadronizerFilter;
1077 DEFINE_FWK_MODULE(Pythia8HepMC3HadronizerFilter);
1078 
1079 typedef edm::ConcurrentGeneratorFilter<Pythia8HepMC3Hadronizer, ConcurrentExternalDecayDriver>
1080     Pythia8HepMC3ConcurrentGeneratorFilter;
1081 DEFINE_FWK_MODULE(Pythia8HepMC3ConcurrentGeneratorFilter);
1082 
1083 typedef edm::ConcurrentHadronizerFilter<Pythia8HepMC3Hadronizer, ConcurrentExternalDecayDriver>
1084     Pythia8HepMC3ConcurrentHadronizerFilter;
1085 DEFINE_FWK_MODULE(Pythia8HepMC3ConcurrentHadronizerFilter);