File indexing completed on 2024-04-25 02:14:04
0001 import FWCore.ParameterSet.Config as cms
0002 from PhysicsTools.NanoAOD.common_cff import *
0003 from PhysicsTools.NanoAOD.simpleGenParticleFlatTableProducer_cfi import simpleGenParticleFlatTableProducer
0004
0005
0006
0007
0008 finalGenParticles = cms.EDProducer("GenParticlePruner",
0009 src = cms.InputTag("prunedGenParticles"),
0010 select = cms.vstring(
0011 "drop *",
0012 "keep++ abs(pdgId) == 15 & (pt > 15 || isPromptDecayed() )",
0013
0014 "keep+ abs(pdgId) == 15 ",
0015 "+keep pdgId == 22 && status == 1 && (pt > 10 || isPromptFinalState())",
0016 "+keep abs(pdgId) == 11 || abs(pdgId) == 13 || abs(pdgId) == 15",
0017 "drop abs(pdgId)= 2212 && abs(pz) > 1000",
0018 "keep (400 < abs(pdgId) < 600) || (4000 < abs(pdgId) < 6000)",
0019 "keep abs(pdgId) == 12 || abs(pdgId) == 14 || abs(pdgId) == 16",
0020 "keep status == 3 || (status > 20 && status < 30)",
0021 "keep isHardProcess() || fromHardProcessDecayed() || fromHardProcessFinalState() || (statusFlags().fromHardProcess() && statusFlags().isLastCopy())",
0022 "keep (status > 70 && status < 80 && pt > 15) ",
0023 "keep abs(pdgId) == 23 || abs(pdgId) == 24 || abs(pdgId) == 25 || abs(pdgId) == 37 ",
0024
0025 "keep (1000001 <= abs(pdgId) <= 1000039 ) || ( 2000001 <= abs(pdgId) <= 2000015)",
0026 )
0027 )
0028
0029
0030
0031
0032 genParticleTable = simpleGenParticleFlatTableProducer.clone(
0033 src = cms.InputTag("finalGenParticles"),
0034 name= cms.string("GenPart"),
0035 doc = cms.string("interesting gen particles "),
0036 externalVariables = cms.PSet(
0037 iso = ExtVar(cms.InputTag("genIso"), float, precision=8, doc="Isolation for leptons"),
0038 ),
0039 variables = cms.PSet(
0040 pt = Var("pt", float, precision=8),
0041 phi = Var("phi", float,precision=8),
0042 eta = Var("eta", float,precision=8),
0043 mass = Var("?!((abs(pdgId)>=1 && abs(pdgId)<=5) || (abs(pdgId)>=11 && abs(pdgId)<=16) || pdgId==21 || pdgId==111 || abs(pdgId)==211 || abs(pdgId)==421 || abs(pdgId)==411 || (pdgId==22 && mass<1))?mass:0", float,precision="?((abs(pdgId)==6 || abs(pdgId)>1000000) && statusFlags().isLastCopy())?20:8",doc="Mass stored for all particles with the exception of quarks (except top), leptons/neutrinos, photons with mass < 1 GeV, gluons, pi0(111), pi+(211), D0(421), and D+(411). For these particles, you can lookup the value from PDG."),
0044 pdgId = Var("pdgId", int, doc="PDG id"),
0045 status = Var("status", int, doc="Particle status. 1=stable"),
0046 genPartIdxMother = Var("?numberOfMothers>0?motherRef(0).key():-1", "int16", doc="index of the mother particle"),
0047 statusFlags = (Var(
0048 "statusFlags().isLastCopyBeforeFSR() * 16384 +"
0049 "statusFlags().isLastCopy() * 8192 +"
0050 "statusFlags().isFirstCopy() * 4096 +"
0051 "statusFlags().fromHardProcessBeforeFSR() * 2048 +"
0052 "statusFlags().isDirectHardProcessTauDecayProduct() * 1024 +"
0053 "statusFlags().isHardProcessTauDecayProduct() * 512 +"
0054 "statusFlags().fromHardProcess() * 256 +"
0055 "statusFlags().isHardProcess() * 128 +"
0056 "statusFlags().isDirectHadronDecayProduct() * 64 +"
0057 "statusFlags().isDirectPromptTauDecayProduct() * 32 +"
0058 "statusFlags().isDirectTauDecayProduct() * 16 +"
0059 "statusFlags().isPromptTauDecayProduct() * 8 +"
0060 "statusFlags().isTauDecayProduct() * 4 +"
0061 "statusFlags().isDecayedLeptonHadron() * 2 +"
0062 "statusFlags().isPrompt() * 1 ",
0063 "uint16", doc=("gen status flags stored bitwise, bits are: "
0064 "0 : isPrompt, "
0065 "1 : isDecayedLeptonHadron, "
0066 "2 : isTauDecayProduct, "
0067 "3 : isPromptTauDecayProduct, "
0068 "4 : isDirectTauDecayProduct, "
0069 "5 : isDirectPromptTauDecayProduct, "
0070 "6 : isDirectHadronDecayProduct, "
0071 "7 : isHardProcess, "
0072 "8 : fromHardProcess, "
0073 "9 : isHardProcessTauDecayProduct, "
0074 "10 : isDirectHardProcessTauDecayProduct, "
0075 "11 : fromHardProcessBeforeFSR, "
0076 "12 : isFirstCopy, "
0077 "13 : isLastCopy, "
0078 "14 : isLastCopyBeforeFSR, ")
0079 )),
0080 )
0081 )
0082
0083 genIso = cms.EDProducer("GenPartIsoProducer",
0084 genPart=cms.InputTag("finalGenParticles"),
0085 packedGenPart=cms.InputTag("packedGenParticles"),
0086 additionalPdgId=cms.int32(22),
0087 )
0088
0089 genParticleTask = cms.Task(finalGenParticles, genIso)
0090 genParticleTablesTask = cms.Task(genParticleTable)