From fe4e75a756afc513dcabf68218229bbd0c4a0721 Mon Sep 17 00:00:00 2001 From: vslokenb <108902248+vslokenb@users.noreply.github.com> Date: Sun, 23 Jul 2023 15:31:53 +0200 Subject: [PATCH 1/6] Add files via upload Plugins needed for recoil to top option. --- .../plugins/Pythia8Hadronizer.cc | 19 +++ .../Pythia8Interface/plugins/TopRecoilHook.h | 136 ++++++++++++++++++ 2 files changed, 155 insertions(+) create mode 100644 GeneratorInterface/Pythia8Interface/plugins/TopRecoilHook.h diff --git a/GeneratorInterface/Pythia8Interface/plugins/Pythia8Hadronizer.cc b/GeneratorInterface/Pythia8Interface/plugins/Pythia8Hadronizer.cc index fbeca29a6c45e..18181e47303e7 100644 --- a/GeneratorInterface/Pythia8Interface/plugins/Pythia8Hadronizer.cc +++ b/GeneratorInterface/Pythia8Interface/plugins/Pythia8Hadronizer.cc @@ -19,6 +19,7 @@ using namespace Pythia8; #include "GeneratorInterface/Pythia8Interface/plugins/ReweightUserHooks.h" #include "GeneratorInterface/Pythia8Interface/interface/CustomHook.h" +#include "GeneratorInterface/Pythia8Interface/plugins/TopRecoilHook.h" // PS matchning prototype // @@ -152,6 +153,9 @@ class Pythia8Hadronizer : public Py8InterfaceBase { //Generic customized hooks vector std::shared_ptr fCustomHooksVector; + //RecoilToTop userhook + std::shared_ptr fTopRecoilHook; + int EV1_nFinal; bool EV1_vetoOn; int EV1_maxVetoCount; @@ -425,6 +429,14 @@ bool Pythia8Hadronizer::initializeForInternalPartons() { (fUserHooksVector->hooks).push_back(fPowhegHooksBB4L); } + bool TopRecoilHook1 = fMasterGen->settings.flag("TopRecoilHook:doTopRecoilIn"); + if (TopRecoilHook1) { + edm::LogInfo("Pythia8Interface") << "Turning on RecoilToTop hook from Pythia8Interface"; + if (!fTopRecoilHook.get()) + fTopRecoilHook.reset(new TopRecoilHook()); + (fUserHooksVector->hooks).push_back(fTopRecoilHook); + } + //adapted from main89.cc in pythia8 examples bool internalMatching = fMasterGen->settings.flag("JetMatching:merge"); bool internalMerging = !(fMasterGen->settings.word("Merging:Process") == "void"); @@ -594,6 +606,13 @@ bool Pythia8Hadronizer::initializeForExternalPartons() { (fUserHooksVector->hooks).push_back(fPowhegHooksBB4L); } + bool TopRecoilHook1 = fMasterGen->settings.flag("TopRecoilHook:doTopRecoilIn"); + if (TopRecoilHook1) { + edm::LogInfo("Pythia8Interface") << "Turning on RecoilToTop hook from Pythia8Interface"; + if (!fTopRecoilHook.get()) + fTopRecoilHook.reset(new TopRecoilHook()); + (fUserHooksVector->hooks).push_back(fTopRecoilHook); + } //adapted from main89.cc in pythia8 examples bool internalMatching = fMasterGen->settings.flag("JetMatching:merge"); bool internalMerging = !(fMasterGen->settings.word("Merging:Process") == "void"); diff --git a/GeneratorInterface/Pythia8Interface/plugins/TopRecoilHook.h b/GeneratorInterface/Pythia8Interface/plugins/TopRecoilHook.h new file mode 100644 index 0000000000000..c07c84ef9e4c9 --- /dev/null +++ b/GeneratorInterface/Pythia8Interface/plugins/TopRecoilHook.h @@ -0,0 +1,136 @@ +// TopRecoilHook.h is a part of the PYTHIA event generator. +// Copyright (C) 2020 Torbjorn Sjostrand. +// PYTHIA is licenced under the GNU GPL v2 or later, see COPYING for details. +// Please respect the MCnet Guidelines, see GUIDELINES for details. + +// Includes a user hook that corrects emission in top decay for dipole +// from gluon to W, to instead be from gluon to top. + +// Important: the top mass shift analysis encoded here is very primitive, +// does not perform well at all, and should not be taken seriously. +// The important part is that you see how the different scenarios +// should be set up to operate as intended. + +#include "Pythia8/Pythia.h" +namespace Pythia8 { + + //========================================================================== + + // Write own derived UserHooks class for modified emission in top decay. + + class TopRecoilHook : public UserHooks { + public: + // Constructor. + // doTopRecoil : eikonal correction in GW dipole on/off when no MEC applied. + // useOldDipole : in GW dipole, use partons before or after branching. + // doList : diagnostic output; set false for production runs. + TopRecoilHook(bool doTopRecoilIn = true, bool useOldDipoleIn = false, bool doListIn = false) { + doTopRecoil = doTopRecoilIn; + useOldDipole = useOldDipoleIn; + doList = doListIn; + // Constructor also creates some histograms for analysis inside User Hook. + wtCorr = new Hist("corrective weight", 100, 0., 2.); + } + + // Destructor prints histogram. + ~TopRecoilHook() override { + if (doTopRecoil) + ; + delete wtCorr; + } + + // Initialise. Only use hook for simple showers with recoilToColoured = off. + bool initAfterBeams() override { + // Switch off if recoilToColoured = on. + bool recoilToColoured = settingsPtr->flag("TimeShower:recoilToColoured"); + if (recoilToColoured) + doTopRecoil = false; + // Flag if W mass term is already accounted for (true) or not (false). + recoilDeadCone = settingsPtr->flag("TimeShower:recoilDeadCone"); + // All ok. + return true; + } + + // Allow a veto after an FSR emission + bool canVetoFSREmission() override { return doTopRecoil; } + + // Access the event after an FSR emission, specifically inside top decay. + bool doVetoFSREmission(int sizeOld, const Event& event, int iSys, bool inResonance) override { + // Check that we are inside a resonance decay. + if (!inResonance) + return false; + + // Check that it is a top decay. + int iTop = partonSystemsPtr->getInRes(iSys); + if (iTop == 0 || event[iTop].idAbs() != 6) + return false; + + // Skip first emission, where ME corrections are already made. + int sizeOut = partonSystemsPtr->sizeOut(iSys); + if (sizeOut == 2) + return false; + + // Location of trial new particles: radiator, emitted, recoiler. + int iRad = sizeOld; + int iEmt = sizeOld + 1; + int iRec = sizeOld + 2; + + // The above partons are after emission; + // alternatively use the ones before. + if (useOldDipole) { + iRad = event[iRad].mother1(); + iRec = event[iRec].mother1(); + } + + // Check if newly emitted gluon matches (anti)top colour line. + if (event[iEmt].id() != 21) + return false; + if (event[iTop].id() == 6) { + if (event[iEmt].col() != event[iTop].col()) + return false; + } else { + if (event[iEmt].acol() != event[iTop].acol()) + return false; + } + + // Recoiler should now be a W, else something is wrong. + if (event[iRec].idAbs() != 24) { + return false; + } + + // Denominator: eikonal weight with W as recoiler. + double pRadRec = event[iRad].p() * event[iRec].p(); + double pRadEmt = event[iRad].p() * event[iEmt].p(); + double pRecEmt = event[iRec].p() * event[iEmt].p(); + double wtW = 2. * pRadRec / (pRadEmt * pRecEmt) - pow2(event[iRad].m() / pRadEmt); + // If recoilDeadCone = on, include W mass term in denominator. + if (recoilDeadCone) + wtW -= pow2(event[iRec].m() / pRecEmt); + + // Numerator: eikonal weight with top as recoiler. + double pRadTop = event[iRad].p() * event[iTop].p(); + double pTopEmt = event[iTop].p() * event[iEmt].p(); + double wtT = + 2. * pRadTop / (pRadEmt * pTopEmt) - pow2(event[iRad].m() / pRadEmt) - pow2(event[iTop].m() / pTopEmt); + + // Histogram weight ratio. + wtCorr->fill(wtT / wtW); + + // List relevant properties. + if (doList) { + partonSystemsPtr->list(); + event.list(); + } + //std::cout << "BEGIN TOP TEST: " << (wtT / wtW) << " | END TOP TEST" << endl; + + // Accept/reject emission. Smooth suppression or step function. + return (wtT < wtW * rndmPtr->flat()); + } + + private: + // Options and Histograms. + bool doTopRecoil, useOldDipole, doList, recoilDeadCone; + Hist* wtCorr; + }; + +} // namespace Pythia8 From 0e6498167ae48f58204d4dc1d037ec190799afbb Mon Sep 17 00:00:00 2001 From: vslokenb <108902248+vslokenb@users.noreply.github.com> Date: Sun, 23 Jul 2023 15:32:37 +0200 Subject: [PATCH 2/6] Add files via upload --- .../Pythia8Interface/src/Py8InterfaceBase.cc | 15 +++++++++++++-- 1 file changed, 13 insertions(+), 2 deletions(-) diff --git a/GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc b/GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc index da9708c4fb253..99f8ec98dfcd1 100644 --- a/GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc +++ b/GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc @@ -70,9 +70,15 @@ namespace gen { } bool Py8InterfaceBase::readSettings(int) { + //Pythia 8's default value for first argument to constructor + const string xmlDir = "../share/Pythia8/xmldoc"; + bool printBanner = true; + if (fParameters.exists("printBanner")) { + printBanner = fParameters.getUntrackedParameter("printBanner"); + } if (!fMasterGen.get()) - fMasterGen = std::make_unique(); - fDecayer = std::make_unique(); + fMasterGen = std::make_unique(xmlDir, printBanner); + fDecayer = std::make_unique(xmlDir, printBanner); //add settings for resonance decay filter fMasterGen->settings.addFlag("BiasedTauDecayer:filter", false); @@ -98,6 +104,11 @@ namespace gen { fMasterGen->settings.addParm("PTFilter:quarkRapidity", 10.0, true, true, 0.0, 10.); fMasterGen->settings.addParm("PTFilter:quarkPt", -.1, true, true, -.1, 100.); + //add settings for RecoilToTop tool + fMasterGen->settings.addFlag("TopRecoilHook:doTopRecoilIn", false); + fMasterGen->settings.addFlag("TopRecoilHook:useOldDipoleIn", false); + fMasterGen->settings.addFlag("TopRecoilHook:doListIn", false); + //add settings for powheg resonance scale calculation fMasterGen->settings.addFlag("POWHEGres:calcScales", false); fMasterGen->settings.addFlag("POWHEG:bb4l", false); From 2cf5c2ff1a2b6eb4f7cf90352058b264b7e1a07b Mon Sep 17 00:00:00 2001 From: vslokenb <108902248+vslokenb@users.noreply.github.com> Date: Sun, 23 Jul 2023 15:33:06 +0200 Subject: [PATCH 3/6] Add files via upload --- .../test/TopRecoilHook_test_cfg.py | 92 +++++++++++++++++++ 1 file changed, 92 insertions(+) create mode 100644 GeneratorInterface/Pythia8Interface/test/TopRecoilHook_test_cfg.py diff --git a/GeneratorInterface/Pythia8Interface/test/TopRecoilHook_test_cfg.py b/GeneratorInterface/Pythia8Interface/test/TopRecoilHook_test_cfg.py new file mode 100644 index 0000000000000..7f8b07665c20c --- /dev/null +++ b/GeneratorInterface/Pythia8Interface/test/TopRecoilHook_test_cfg.py @@ -0,0 +1,92 @@ +import sys +import FWCore.ParameterSet.Config as cms +import FWCore.ParameterSet.VarParsing as VarParsing +#import pythia8 + +options = VarParsing.VarParsing ('standard') +options.register('runOnly', '', VarParsing.VarParsing.multiplicity.singleton,VarParsing.VarParsing.varType.string, "Run only specified analysis") +options.register('yodafile', 'test_toprecoil.yoda', VarParsing.VarParsing.multiplicity.singleton,VarParsing.VarParsing.varType.string, "Name of yoda output file") +options.register('recoiltoTop','on',VarParsing.VarParsing.multiplicity.singleton,VarParsing.VarParsing.varType.string, "TopRecoilHook setting") +options.register('recoiltoBottom','off',VarParsing.VarParsing.multiplicity.singleton,VarParsing.VarParsing.varType.string, "RecoilToColor setting") +options.setDefault('maxEvents', 1000) +options.register('alphaSvalue', 0.1365,VarParsing.VarParsing.multiplicity.singleton,VarParsing.VarParsing.varType.string, "AlphaS value setting") +options.register('topmass', 172.5,VarParsing.VarParsing.multiplicity.singleton,VarParsing.VarParsing.varType.string, "Top mass setting") +if(hasattr(sys, "argv")): + options.parseArguments() +print(options) + +process = cms.Process("runRivetAnalysis") + +from Configuration.Generator.Pythia8CommonSettings_cfi import * +from Configuration.Generator.Pythia8CUEP8M1Settings_cfi import * +from Configuration.Generator.Pythia8aMCatNLOSettings_cfi import * +from IOMC.RandomEngine.RandomServiceHelper import RandomNumberServiceHelper +#randSvc = RandomNumberServiceHelper(process.RandomNumberGeneratorService) +#randSvc.populate() + +# import of standard configurations +process.load('Configuration.StandardSequences.Services_cff') + +randSvc = RandomNumberServiceHelper(process.RandomNumberGeneratorService) +randSvc.populate() + +process.load("FWCore.MessageLogger.MessageLogger_cfi") +process.MessageLogger.cerr.FwkReport.reportEvery = cms.untracked.int32(1000) +process.maxEvents = cms.untracked.PSet(input = cms.untracked.int32(options.maxEvents)) + +process.source = cms.Source("EmptySource") + +process.load('SimGeneral.HepPDTESSource.pythiapdt_cfi') +process.generator = cms.EDFilter("Pythia8GeneratorFilter", + comEnergy = cms.double(13000.0), + crossSection = cms.untracked.double(421.1), + filterEfficiency = cms.untracked.double(1), + maxEventsToPrint = cms.untracked.int32(0), + pythiaHepMCVerbosity = cms.untracked.bool(False), + pythiaPylistVerbosity = cms.untracked.int32(1), + PythiaParameters = cms.PSet( + pythia8CommonSettingsBlock, + pythia8CUEP8M1SettingsBlock, + pythia8aMCatNLOSettingsBlock, + processParameters = cms.vstring( + 'Main:timesAllowErrors = 10000', + 'ParticleDecays:limitTau0 = on', + 'ParticleDecays:tauMax = 10', + 'Tune:ee 7', + 'Tune:pp 14', # Monash tune + 'Top:gg2ttbar = on', + 'Top:qqbar2ttbar = on', + '6:m0 = %s'%options.topmass, # top mass' + 'TopRecoilHook:doTopRecoilIn = %s'%options.recoiltoTop, + 'TimeShower:recoilToColoured = %s'%options.recoiltoBottom, + 'TimeShower:alphaSvalue = %s' %options.alphaSvalue + ), + parameterSets = cms.vstring('processParameters',) + ) +) + +process.load("GeneratorInterface.RivetInterface.rivetAnalyzer_cfi") + +#pythia.Settings.listAll() + +if options.runOnly: + process.rivetAnalyzer.AnalysisNames = cms.vstring(options.runOnly) +else: + process.rivetAnalyzer.AnalysisNames = cms.vstring( + 'CMS_2018_I1663958', # diff xs lepton+jets + 'CMS_2018_I1662081', # event variables lepton+jets + 'MC_TOPMASS_LJETS', # MC analysis for lepton+jets top mass + 'CMS_2018_I1690148', # jet substructure + 'MC_BFRAG_LJETS', # b fragmentation + 'CMS_2018_I1703993', # diff xs dilepton + 'CMS_2019_I1764472', #boosted top mass + 'MC_TOPMASS_LSMT', #dilepton invariant mass, soft muon +primary + 'ATLAS_2019_I1724098', #missing momentum analysis + ) +process.rivetAnalyzer.OutputFile = options.yodafile +process.rivetAnalyzer.HepMCCollection = cms.InputTag("generator:unsmeared") +process.rivetAnalyzer.CrossSection = 831.76 # NNLO (arXiv:1303.6254) + +process.p = cms.Path(process.generator*process.rivetAnalyzer) + + From 4d0724fb9528d8f9e41f74b092ad6311291ea17a Mon Sep 17 00:00:00 2001 From: vslokenb <108902248+vslokenb@users.noreply.github.com> Date: Sun, 23 Jul 2023 15:37:44 +0200 Subject: [PATCH 4/6] Add files via upload --- GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc b/GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc index 99f8ec98dfcd1..d0cd776eaeb6d 100644 --- a/GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc +++ b/GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc @@ -77,8 +77,8 @@ namespace gen { printBanner = fParameters.getUntrackedParameter("printBanner"); } if (!fMasterGen.get()) - fMasterGen = std::make_unique(xmlDir, printBanner); - fDecayer = std::make_unique(xmlDir, printBanner); + fMasterGen = std::make_unique(); + fDecayer = std::make_unique(); //add settings for resonance decay filter fMasterGen->settings.addFlag("BiasedTauDecayer:filter", false); From bea808a0df1e17a5aa42f3863bb80f0218961755 Mon Sep 17 00:00:00 2001 From: vslokenb <108902248+vslokenb@users.noreply.github.com> Date: Sun, 23 Jul 2023 15:40:20 +0200 Subject: [PATCH 5/6] Add files via upload --- .../Pythia8Interface/src/Py8InterfaceBase.cc | 8 +------- 1 file changed, 1 insertion(+), 7 deletions(-) diff --git a/GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc b/GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc index d0cd776eaeb6d..1a524e794d83a 100644 --- a/GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc +++ b/GeneratorInterface/Pythia8Interface/src/Py8InterfaceBase.cc @@ -69,13 +69,7 @@ namespace gen { } } - bool Py8InterfaceBase::readSettings(int) { - //Pythia 8's default value for first argument to constructor - const string xmlDir = "../share/Pythia8/xmldoc"; - bool printBanner = true; - if (fParameters.exists("printBanner")) { - printBanner = fParameters.getUntrackedParameter("printBanner"); - } + bool Py8InterfaceBase::readSettings(int) { if (!fMasterGen.get()) fMasterGen = std::make_unique(); fDecayer = std::make_unique(); From 5a9ed851f6717ca84120e81106ae27755882c301 Mon Sep 17 00:00:00 2001 From: vslokenb <108902248+vslokenb@users.noreply.github.com> Date: Thu, 27 Jul 2023 15:24:08 +0200 Subject: [PATCH 6/6] Update TopRecoilHook.h Remove useless check as suggested by PR 42313. --- GeneratorInterface/Pythia8Interface/plugins/TopRecoilHook.h | 6 +----- 1 file changed, 1 insertion(+), 5 deletions(-) diff --git a/GeneratorInterface/Pythia8Interface/plugins/TopRecoilHook.h b/GeneratorInterface/Pythia8Interface/plugins/TopRecoilHook.h index c07c84ef9e4c9..309089939325c 100644 --- a/GeneratorInterface/Pythia8Interface/plugins/TopRecoilHook.h +++ b/GeneratorInterface/Pythia8Interface/plugins/TopRecoilHook.h @@ -33,11 +33,7 @@ namespace Pythia8 { } // Destructor prints histogram. - ~TopRecoilHook() override { - if (doTopRecoil) - ; - delete wtCorr; - } + ~TopRecoilHook() override { delete wtCorr; } // Initialise. Only use hook for simple showers with recoilToColoured = off. bool initAfterBeams() override {