diff --git a/Configuration/Generator/python/Py8Eta2MuGammaPtExpGun_cfi.py b/Configuration/Generator/python/Py8Eta2MuGammaPtExpGun_cfi.py new file mode 100644 index 0000000000000..03febe3e40912 --- /dev/null +++ b/Configuration/Generator/python/Py8Eta2MuGammaPtExpGun_cfi.py @@ -0,0 +1,47 @@ +import FWCore.ParameterSet.Config as cms + +from Configuration.Generator.Pythia8CommonSettings_cfi import * +from Configuration.Generator.MCTunes2017.PythiaCP5Settings_cfi import * + +generator = cms.EDFilter("Pythia8PtExpGun", + + maxEventsToPrint = cms.untracked.int32(1), + pythiaPylistVerbosity = cms.untracked.int32(1), + pythiaHepMCVerbosity = cms.untracked.bool(True), + + PGunParameters = cms.PSet( + ParticleID = cms.vint32(221), + AddAntiParticle = cms.bool(False), + MinPhi = cms.double(-3.14159265359), + MaxPhi = cms.double(3.14159265359), + #MinPt = cms.double(5.0), + MinPt = cms.double(10.0), + #MaxPt = cms.double(65.0), + #MaxPt = cms.double(25.0), + MaxPt = cms.double(65.0), + MinEta = cms.double(-2.4), + MaxEta = cms.double(2.4) + ), + + PythiaParameters = cms.PSet( + pythia8CommonSettingsBlock, + pythia8CP5SettingsBlock, + processParameters = cms.vstring( + #'SLHA:keepSM = on', + #'SLHA:minMassSM = 10.', + # Very important to enable override! + 'SLHA:allowUserOverride = on', + 'RHadrons:allow = on', + 'RHadrons:allowDecay = on', + #'32:mayDecay = true', + '221:mayDecay = true', + # Set decay channels of eta (mumugamma) + '221:oneChannel = 1 1.0 0 13 -13 22' + ), + parameterSets = cms.vstring( + 'pythia8CommonSettings', + 'pythia8CP5Settings', + 'processParameters', + ) + ) +) diff --git a/GeneratorInterface/Pythia8Interface/plugins/Py8PtExpGun.cc b/GeneratorInterface/Pythia8Interface/plugins/Py8PtExpGun.cc new file mode 100644 index 0000000000000..9ac3ff1cc2778 --- /dev/null +++ b/GeneratorInterface/Pythia8Interface/plugins/Py8PtExpGun.cc @@ -0,0 +1,124 @@ +#include + +#include "GeneratorInterface/Core/interface/GeneratorFilter.h" +#include "GeneratorInterface/ExternalDecays/interface/ExternalDecayDriver.h" + +#include "GeneratorInterface/Pythia8Interface/interface/Py8GunBase.h" + +namespace gen { + + class Py8PtExpGun : public Py8GunBase { + public: + Py8PtExpGun(edm::ParameterSet const&); + ~Py8PtExpGun() override {} + + bool generatePartonsAndHadronize() override; + const char* classname() const override; + + private: + // PtExpGun particle(s) characteristics + double fMinEta; + double fMaxEta; + double fMinPt; + double fMaxPt; + bool fAddAntiParticle; + }; + + // implementation + // + Py8PtExpGun::Py8PtExpGun(edm::ParameterSet const& ps) : Py8GunBase(ps) { + // ParameterSet defpset ; + edm::ParameterSet pgun_params = ps.getParameter("PGunParameters"); // , defpset ) ; + fMinEta = pgun_params.getParameter("MinEta"); // ,-2.2); + fMaxEta = pgun_params.getParameter("MaxEta"); // , 2.2); + fMinPt = pgun_params.getParameter("MinPt"); // , 0.); + fMaxPt = pgun_params.getParameter("MaxPt"); // , 0.); + fAddAntiParticle = pgun_params.getParameter("AddAntiParticle"); //, false) ; + } + + bool Py8PtExpGun::generatePartonsAndHadronize() { + fMasterGen->event.reset(); + + int NTotParticles = fPartIDs.size(); + if (fAddAntiParticle) + NTotParticles *= 2; + + // energy below is dummy, it is not used + (fMasterGen->event).append(990, -11, 0, 0, 2, 1 + NTotParticles, 0, 0, 0., 0., 0., 15000., 15000.); + + for (size_t i = 0; i < fPartIDs.size(); i++) { + int particleID = fPartIDs[i]; // this is PDG - need to convert to Py8 ??? + + double phi = (fMaxPhi - fMinPhi) * randomEngine().flat() + fMinPhi; + double eta = (fMaxEta - fMinEta) * randomEngine().flat() + fMinEta; + double the = 2. * atan(exp(-eta)); + + //-log(flat) = exponential distribution + // need the /10.0 and the min with 1.0 to make sure pt doesn't go too high + // 10.0 chosen to give last pt bin (overflow bin) a reasonable (not unnaturally high) content + double pt = (std::min(-1 / 10.0 * log(randomEngine().flat()), 1.0)) * (fMaxPt - fMinPt) + fMinPt; + + double mass = (fMasterGen->particleData).m0(particleID); + + double pp = pt / sin(the); // sqrt( ee*ee - mass*mass ); + double ee = sqrt(pp * pp + mass * mass); + + double px = pt * cos(phi); + double py = pt * sin(phi); + double pz = pp * cos(the); + + if (!((fMasterGen->particleData).isParticle(particleID))) { + particleID = std::abs(particleID); + } + if (1 <= std::abs(particleID) && std::abs(particleID) <= 6) // quarks + (fMasterGen->event).append(particleID, 23, 1, 0, 0, 0, 101, 0, px, py, pz, ee, mass); + else if (std::abs(particleID) == 21) // gluons + (fMasterGen->event).append(21, 23, 1, 0, 0, 0, 101, 102, px, py, pz, ee, mass); + // other + else { + (fMasterGen->event).append(particleID, 1, 1, 0, 0, 0, 0, 0, px, py, pz, ee, mass); + int eventSize = (fMasterGen->event).size() - 1; + // -log(flat) = exponential distribution + double tauTmp = -(fMasterGen->event)[eventSize].tau0() * log(randomEngine().flat()); + (fMasterGen->event)[eventSize].tau(tauTmp); + } + + // Here also need to add anti-particle (if any) + // otherwise just add a 2nd particle of the same type + // (for example, gamma) + // + if (fAddAntiParticle) { + if (1 <= std::abs(particleID) && std::abs(particleID) <= 6) { // quarks + (fMasterGen->event).append(-particleID, 23, 1, 0, 0, 0, 0, 101, -px, -py, -pz, ee, mass); + } else if (std::abs(particleID) == 21) { // gluons + (fMasterGen->event).append(21, 23, 1, 0, 0, 0, 102, 101, -px, -py, -pz, ee, mass); + } else { + if ((fMasterGen->particleData).isParticle(-particleID)) { + (fMasterGen->event).append(-particleID, 1, 1, 0, 0, 0, 0, 0, -px, -py, -pz, ee, mass); + } else { + (fMasterGen->event).append(particleID, 1, 1, 0, 0, 0, 0, 0, -px, -py, -pz, ee, mass); + } + int eventSize = (fMasterGen->event).size() - 1; + // -log(flat) = exponential distribution + double tauTmp = -(fMasterGen->event)[eventSize].tau0() * log(randomEngine().flat()); + (fMasterGen->event)[eventSize].tau(tauTmp); + } + } + } + + if (!fMasterGen->next()) + return false; + evtGenDecay(); + + event() = std::make_unique(); + return toHepMC.fill_next_event(fMasterGen->event, event().get()); + } + + const char* Py8PtExpGun::classname() const { return "Py8PtExpGun"; } + + typedef edm::GeneratorFilter Pythia8PtExpGun; + +} // namespace gen + +using gen::Pythia8PtExpGun; +DEFINE_FWK_MODULE(Pythia8PtExpGun);