CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
jets_cff.py
Go to the documentation of this file.
1 import FWCore.ParameterSet.Config as cms
2 from Configuration.Eras.Modifier_run2_jme_2016_cff import run2_jme_2016
3 from Configuration.Eras.Modifier_run2_jme_2017_cff import run2_jme_2017
4 from Configuration.Eras.Modifier_run2_miniAOD_80XLegacy_cff import run2_miniAOD_80XLegacy
5 
8 
9 
10 ##################### User floats producers, selectors ##########################
11 from RecoJets.JetProducers.ak4PFJets_cfi import ak4PFJets
12 
13 chsForSATkJets = cms.EDFilter("CandPtrSelector", src = cms.InputTag("packedPFCandidates"), cut = cms.string('charge()!=0 && pvAssociationQuality()>=5 && vertexRef().key()==0'))
14 softActivityJets = ak4PFJets.clone(src = 'chsForSATkJets', doAreaFastjet = False, jetPtMin=1)
15 softActivityJets10 = cms.EDFilter("CandPtrSelector", src = cms.InputTag("softActivityJets"), cut = cms.string('pt>10'))
16 softActivityJets5 = cms.EDFilter("CandPtrSelector", src = cms.InputTag("softActivityJets"), cut = cms.string('pt>5'))
17 softActivityJets2 = cms.EDFilter("CandPtrSelector", src = cms.InputTag("softActivityJets"), cut = cms.string('pt>2'))
18 
20 # Note: Safe to always add 'L2L3Residual' as MC contains dummy L2L3Residual corrections (always set to 1)
21 # (cf. https://twiki.cern.ch/twiki/bin/view/CMSPublic/WorkBookJetEnergyCorrections#CMSSW_7_6_4_and_above )
22 jetCorrFactorsNano = patJetCorrFactors.clone(src='slimmedJets',
23  levels = cms.vstring('L1FastJet',
24  'L2Relative',
25  'L3Absolute',
26  'L2L3Residual'),
27  primaryVertices = cms.InputTag("offlineSlimmedPrimaryVertices"),
28 )
29 jetCorrFactorsAK8 = patJetCorrFactors.clone(src='slimmedJetsAK8',
30  levels = cms.vstring('L1FastJet',
31  'L2Relative',
32  'L3Absolute',
33  'L2L3Residual'),
34  payload = cms.string('AK8PFPuppi'),
35  primaryVertices = cms.InputTag("offlineSlimmedPrimaryVertices"),
36 )
37 run2_miniAOD_80XLegacy.toModify(jetCorrFactorsAK8, payload = cms.string('AK8PFchs')) # ak8PFJetsCHS in 2016 80X miniAOD
38 
40 
41 updatedJets = updatedPatJets.clone(
42  addBTagInfo=False,
43  jetSource='slimmedJets',
44  jetCorrFactorsSource=cms.VInputTag(cms.InputTag("jetCorrFactorsNano") ),
45 )
46 
47 updatedJetsAK8 = updatedPatJets.clone(
48  addBTagInfo=False,
49  jetSource='slimmedJetsAK8',
50  jetCorrFactorsSource=cms.VInputTag(cms.InputTag("jetCorrFactorsAK8") ),
51 )
52 
53 
54 looseJetId = cms.EDProducer("PatJetIDValueMapProducer",
55  filterParams=cms.PSet(
56  version = cms.string('WINTER16'),
57  quality = cms.string('LOOSE'),
58  ),
59  src = cms.InputTag("updatedJets")
60 )
61 tightJetId = cms.EDProducer("PatJetIDValueMapProducer",
62  filterParams=cms.PSet(
63  version = cms.string('SUMMER18'),
64  quality = cms.string('TIGHT'),
65  ),
66  src = cms.InputTag("updatedJets")
67 )
68 tightJetIdLepVeto = cms.EDProducer("PatJetIDValueMapProducer",
69  filterParams=cms.PSet(
70  version = cms.string('SUMMER18'),
71  quality = cms.string('TIGHTLEPVETO'),
72  ),
73  src = cms.InputTag("updatedJets")
74 )
75 run2_jme_2016.toModify( tightJetId.filterParams, version = "WINTER16" )
76 run2_jme_2016.toModify( tightJetIdLepVeto.filterParams, version = "WINTER16" )
77 run2_jme_2017.toModify( tightJetId.filterParams, version = "WINTER17" )
78 run2_jme_2017.toModify( tightJetIdLepVeto.filterParams, version = "WINTER17" )
79 
80 
81 looseJetIdAK8 = cms.EDProducer("PatJetIDValueMapProducer",
82  filterParams=cms.PSet(
83  version = cms.string('WINTER16'),
84  quality = cms.string('LOOSE'),
85  ),
86  src = cms.InputTag("updatedJetsAK8")
87 )
88 tightJetIdAK8 = cms.EDProducer("PatJetIDValueMapProducer",
89  filterParams=cms.PSet(
90  version = cms.string('SUMMER18PUPPI'),
91  quality = cms.string('TIGHT'),
92  ),
93  src = cms.InputTag("updatedJetsAK8")
94 )
95 tightJetIdLepVetoAK8 = cms.EDProducer("PatJetIDValueMapProducer",
96  filterParams=cms.PSet(
97  version = cms.string('SUMMER18PUPPI'),
98  quality = cms.string('TIGHTLEPVETO'),
99  ),
100  src = cms.InputTag("updatedJetsAK8")
101 )
102 run2_jme_2016.toModify( tightJetIdAK8.filterParams, version = "WINTER16" )
103 run2_jme_2016.toModify( tightJetIdLepVetoAK8.filterParams, version = "WINTER16" )
104 run2_jme_2017.toModify( tightJetIdAK8.filterParams, version = "WINTER17PUPPI" )
105 run2_jme_2017.toModify( tightJetIdLepVetoAK8.filterParams, version = "WINTER17PUPPI" )
106 
107 
108 bJetVars = cms.EDProducer("JetRegressionVarProducer",
109  pvsrc = cms.InputTag("offlineSlimmedPrimaryVertices"),
110  src = cms.InputTag("updatedJets"),
111  svsrc = cms.InputTag("slimmedSecondaryVertices"),
112  gpsrc = cms.InputTag("prunedGenParticles"),
113  #musrc = cms.InputTag("slimmedMuons"),
114  #elesrc = cms.InputTag("slimmedElectrons")
115 )
116 
117 jercVars = cms.EDProducer("BetaStarPackedCandidateVarProducer",
118  srcJet = cms.InputTag("updatedJets"),
119  srcPF = cms.InputTag("packedPFCandidates"),
120  maxDR = cms.double(0.4)
121 )
122 
123 
124 updatedJetsWithUserData = cms.EDProducer("PATJetUserDataEmbedder",
125  src = cms.InputTag("updatedJets"),
126  userFloats = cms.PSet(
127  leadTrackPt = cms.InputTag("bJetVars:leadTrackPt"),
128  leptonPtRel = cms.InputTag("bJetVars:leptonPtRel"),
129  leptonPtRatio = cms.InputTag("bJetVars:leptonPtRatio"),
130  leptonPtRelInv = cms.InputTag("bJetVars:leptonPtRelInv"),
131  leptonPtRelv0 = cms.InputTag("bJetVars:leptonPtRelv0"),
132  leptonPtRatiov0 = cms.InputTag("bJetVars:leptonPtRatiov0"),
133  leptonPtRelInvv0 = cms.InputTag("bJetVars:leptonPtRelInvv0"),
134  leptonDeltaR = cms.InputTag("bJetVars:leptonDeltaR"),
135  leptonPt = cms.InputTag("bJetVars:leptonPt"),
136  vtxPt = cms.InputTag("bJetVars:vtxPt"),
137  vtxMass = cms.InputTag("bJetVars:vtxMass"),
138  vtx3dL = cms.InputTag("bJetVars:vtx3dL"),
139  vtx3deL = cms.InputTag("bJetVars:vtx3deL"),
140  ptD = cms.InputTag("bJetVars:ptD"),
141  genPtwNu = cms.InputTag("bJetVars:genPtwNu"),
142  qgl = cms.InputTag('qgtagger:qgLikelihood'),
143  jercCHPUF = cms.InputTag("jercVars:chargedHadronPUEnergyFraction"),
144  jercCHF = cms.InputTag("jercVars:chargedHadronCHSEnergyFraction"),
145  ),
146  userInts = cms.PSet(
147  tightId = cms.InputTag("tightJetId"),
148  tightIdLepVeto = cms.InputTag("tightJetIdLepVeto"),
149  vtxNtrk = cms.InputTag("bJetVars:vtxNtrk"),
150  leptonPdgId = cms.InputTag("bJetVars:leptonPdgId"),
151 
152  ),
153 )
154 run2_jme_2016.toModify(updatedJetsWithUserData.userInts,
155  looseId = cms.InputTag("looseJetId"),
156 )
157 
158 updatedJetsAK8WithUserData = cms.EDProducer("PATJetUserDataEmbedder",
159  src = cms.InputTag("updatedJetsAK8"),
160  userFloats = cms.PSet(),
161  userInts = cms.PSet(
162  tightId = cms.InputTag("tightJetIdAK8"),
163  tightIdLepVeto = cms.InputTag("tightJetIdLepVetoAK8"),
164  ),
165 )
166 run2_jme_2016.toModify(updatedJetsAK8WithUserData.userInts,
167  looseId = cms.InputTag("looseJetIdAK8"),
168 )
169 
170 
171 finalJets = cms.EDFilter("PATJetRefSelector",
172  src = cms.InputTag("updatedJetsWithUserData"),
173  cut = cms.string("pt > 15")
174 )
175 
176 finalJetsAK8 = cms.EDFilter("PATJetRefSelector",
177  src = cms.InputTag("updatedJetsAK8WithUserData"),
178  cut = cms.string("pt > 170")
179 )
180 
181 
182 
183 
184 
185 ##################### Tables for final output and docs ##########################
186 
187 
188 
189 jetTable = cms.EDProducer("SimpleCandidateFlatTableProducer",
190  src = cms.InputTag("linkedObjects","jets"),
191  cut = cms.string(""), #we should not filter on cross linked collections
192  name = cms.string("Jet"),
193  doc = cms.string("slimmedJets, i.e. ak4 PFJets CHS with JECs applied, after basic selection (" + finalJets.cut.value()+")"),
194  singleton = cms.bool(False), # the number of entries is variable
195  extension = cms.bool(False), # this is the main table for the jets
196  externalVariables = cms.PSet(
197  #bRegOld = ExtVar(cms.InputTag("bjetMVA"),float, doc="pt corrected with b-jet regression",precision=14),
198  bRegCorr = ExtVar(cms.InputTag("bjetNN:corr"),float, doc="pt correction for b-jet energy regression",precision=12),
199  bRegRes = ExtVar(cms.InputTag("bjetNN:res"),float, doc="res on pt corrected with b-jet regression",precision=8),
200  ),
201  variables = cms.PSet(P4Vars,
202  area = Var("jetArea()", float, doc="jet catchment area, for JECs",precision=10),
203  nMuons = Var("?hasOverlaps('muons')?overlaps('muons').size():0", int, doc="number of muons in the jet"),
204  muonIdx1 = Var("?overlaps('muons').size()>0?overlaps('muons')[0].key():-1", int, doc="index of first matching muon"),
205  muonIdx2 = Var("?overlaps('muons').size()>1?overlaps('muons')[1].key():-1", int, doc="index of second matching muon"),
206  electronIdx1 = Var("?overlaps('electrons').size()>0?overlaps('electrons')[0].key():-1", int, doc="index of first matching electron"),
207  electronIdx2 = Var("?overlaps('electrons').size()>1?overlaps('electrons')[1].key():-1", int, doc="index of second matching electron"),
208  nElectrons = Var("?hasOverlaps('electrons')?overlaps('electrons').size():0", int, doc="number of electrons in the jet"),
209  btagCMVA = Var("bDiscriminator('pfCombinedMVAV2BJetTags')",float,doc="CMVA V2 btag discriminator",precision=10),
210  btagDeepB = Var("bDiscriminator('pfDeepCSVJetTags:probb')+bDiscriminator('pfDeepCSVJetTags:probbb')",float,doc="DeepCSV b+bb tag discriminator",precision=10),
211  btagDeepFlavB = Var("bDiscriminator('pfDeepFlavourJetTags:probb')+bDiscriminator('pfDeepFlavourJetTags:probbb')+bDiscriminator('pfDeepFlavourJetTags:problepb')",float,doc="DeepFlavour b+bb+lepb tag discriminator",precision=10),
212  btagCSVV2 = Var("bDiscriminator('pfCombinedInclusiveSecondaryVertexV2BJetTags')",float,doc=" pfCombinedInclusiveSecondaryVertexV2 b-tag discriminator (aka CSVV2)",precision=10),
213  btagDeepC = Var("bDiscriminator('pfDeepCSVJetTags:probc')",float,doc="DeepCSV charm btag discriminator",precision=10),
214  btagDeepFlavC = Var("bDiscriminator('pfDeepFlavourJetTags:probc')",float,doc="DeepFlavour charm tag discriminator",precision=10),
215  #puIdDisc = Var("userFloat('pileupJetId:fullDiscriminant')",float,doc="Pilup ID discriminant",precision=10),
216  puId = Var("userInt('pileupJetId:fullId')",int,doc="Pilup ID flags"),
217  jetId = Var("userInt('tightId')*2+4*userInt('tightIdLepVeto')",int,doc="Jet ID flags bit1 is loose (always false in 2017 since it does not exist), bit2 is tight, bit3 is tightLepVeto"),
218  qgl = Var("userFloat('qgl')",float,doc="Quark vs Gluon likelihood discriminator",precision=10),
219  nConstituents = Var("numberOfDaughters()",int,doc="Number of particles in the jet"),
220  rawFactor = Var("1.-jecFactor('Uncorrected')",float,doc="1 - Factor to get back to raw pT",precision=6),
221  chHEF = Var("chargedHadronEnergyFraction()", float, doc="charged Hadron Energy Fraction", precision= 6),
222  neHEF = Var("neutralHadronEnergyFraction()", float, doc="neutral Hadron Energy Fraction", precision= 6),
223  chEmEF = Var("chargedEmEnergyFraction()", float, doc="charged Electromagnetic Energy Fraction", precision= 6),
224  neEmEF = Var("neutralEmEnergyFraction()", float, doc="neutral Electromagnetic Energy Fraction", precision= 6),
225  muEF = Var("muonEnergyFraction()", float, doc="muon Energy Fraction", precision= 6),
226  jercCHPUF = Var("userFloat('jercCHPUF')", float, doc="Pileup Charged Hadron Energy Fraction with the JERC group definition", precision= 6),
227  jercCHF = Var("userFloat('jercCHF')", float, doc="Charged Hadron Energy Fraction with the JERC group definition", precision= 6),
228  )
229 )
230 
231 #jets are not as precise as muons
232 jetTable.variables.pt.precision=10
233 
234 ### Era dependent customization
235 run2_jme_2016.toModify( jetTable.variables, jetId = Var("userInt('tightIdLepVeto')*4+userInt('tightId')*2+userInt('looseId')",int,doc="Jet ID flags bit1 is loose, bit2 is tight, bit3 is tightLepVeto"))
236 
237 bjetMVA= cms.EDProducer("BJetEnergyRegressionMVA",
238  backend = cms.string("TMVA"),
239  src = cms.InputTag("linkedObjects","jets"),
240  pvsrc = cms.InputTag("offlineSlimmedPrimaryVertices"),
241  svsrc = cms.InputTag("slimmedSecondaryVertices"),
242  rhosrc = cms.InputTag("fixedGridRhoFastjetAll"),
243  weightFile = cms.FileInPath("PhysicsTools/NanoAOD/data/bjet-regression.xml"),
244  name = cms.string("JetReg"),
245  isClassifier = cms.bool(False),
246  variablesOrder = cms.vstring(["Jet_pt","nPVs","Jet_eta","Jet_mt","Jet_leadTrackPt","Jet_leptonPtRel","Jet_leptonPt","Jet_leptonDeltaR","Jet_neHEF","Jet_neEmEF","Jet_vtxPt","Jet_vtxMass","Jet_vtx3dL","Jet_vtxNtrk","Jet_vtx3deL"]),
247  variables = cms.PSet(
248  Jet_pt = cms.string("pt"),
249  Jet_eta = cms.string("eta"),
250  Jet_mt = cms.string("mt"),
251  Jet_leadTrackPt = cms.string("userFloat('leadTrackPt')"),
252  Jet_vtxNtrk = cms.string("userInt('vtxNtrk')"),
253  Jet_vtxMass = cms.string("userFloat('vtxMass')"),
254  Jet_vtx3dL = cms.string("userFloat('vtx3dL')"),
255  Jet_vtx3deL = cms.string("userFloat('vtx3deL')"),
256  Jet_vtxPt = cms.string("userFloat('vtxPt')"),
257  Jet_leptonPtRel = cms.string("userFloat('leptonPtRelv0')"),
258  Jet_leptonPt = cms.string("?overlaps('muons').size()>0?overlaps('muons')[0].pt():(?overlaps('electrons').size()>0?overlaps('electrons')[0].pt():0)"),
259  Jet_neHEF = cms.string("neutralHadronEnergy()/energy()"),
260  Jet_neEmEF = cms.string("neutralEmEnergy()/energy()"),
261  Jet_leptonDeltaR = cms.string('''?overlaps('muons').size()>0?deltaR(eta,phi,overlaps('muons')[0].eta,overlaps('muons')[0].phi):
262  (?overlaps('electrons').size()>0?deltaR(eta,phi,overlaps('electrons')[0].eta,overlaps('electrons')[0].phi):
263  0)'''),
264  )
265 
266 )
267 
268 bjetNN= cms.EDProducer("BJetEnergyRegressionMVA",
269  backend = cms.string("TF"),
270  src = cms.InputTag("linkedObjects","jets"),
271  pvsrc = cms.InputTag("offlineSlimmedPrimaryVertices"),
272  svsrc = cms.InputTag("slimmedSecondaryVertices"),
273  rhosrc = cms.InputTag("fixedGridRhoFastjetAll"),
274 
275  weightFile = cms.FileInPath("PhysicsTools/NanoAOD/data/breg_training_2017.pb"),
276  name = cms.string("JetRegNN"),
277  isClassifier = cms.bool(False),
278  variablesOrder = cms.vstring(["Jet_pt","Jet_eta","rho","Jet_mt","Jet_leadTrackPt","Jet_leptonPtRel","Jet_leptonDeltaR","Jet_neHEF","Jet_neEmEF","Jet_vtxPt","Jet_vtxMass","Jet_vtx3dL","Jet_vtxNtrk","Jet_vtx3deL","Jet_numDaughters_pt03","Jet_energyRing_dR0_em_Jet_rawEnergy","Jet_energyRing_dR1_em_Jet_rawEnergy","Jet_energyRing_dR2_em_Jet_rawEnergy","Jet_energyRing_dR3_em_Jet_rawEnergy","Jet_energyRing_dR4_em_Jet_rawEnergy","Jet_energyRing_dR0_neut_Jet_rawEnergy","Jet_energyRing_dR1_neut_Jet_rawEnergy","Jet_energyRing_dR2_neut_Jet_rawEnergy","Jet_energyRing_dR3_neut_Jet_rawEnergy","Jet_energyRing_dR4_neut_Jet_rawEnergy","Jet_energyRing_dR0_ch_Jet_rawEnergy","Jet_energyRing_dR1_ch_Jet_rawEnergy","Jet_energyRing_dR2_ch_Jet_rawEnergy","Jet_energyRing_dR3_ch_Jet_rawEnergy","Jet_energyRing_dR4_ch_Jet_rawEnergy","Jet_energyRing_dR0_mu_Jet_rawEnergy","Jet_energyRing_dR1_mu_Jet_rawEnergy","Jet_energyRing_dR2_mu_Jet_rawEnergy","Jet_energyRing_dR3_mu_Jet_rawEnergy","Jet_energyRing_dR4_mu_Jet_rawEnergy","Jet_chHEF","Jet_chEmEF","Jet_leptonPtRelInv","isEle","isMu","isOther","Jet_mass","Jet_ptd"]),
279  variables = cms.PSet(
280  Jet_pt = cms.string("pt*jecFactor('Uncorrected')"),
281  Jet_mt = cms.string("mt*jecFactor('Uncorrected')"),
282  Jet_eta = cms.string("eta"),
283  Jet_mass = cms.string("mass*jecFactor('Uncorrected')"),
284  Jet_ptd = cms.string("userFloat('ptD')"),
285  Jet_leadTrackPt = cms.string("userFloat('leadTrackPt')"),
286  Jet_vtxNtrk = cms.string("userInt('vtxNtrk')"),
287  Jet_vtxMass = cms.string("userFloat('vtxMass')"),
288  Jet_vtx3dL = cms.string("userFloat('vtx3dL')"),
289  Jet_vtx3deL = cms.string("userFloat('vtx3deL')"),
290  Jet_vtxPt = cms.string("userFloat('vtxPt')"),
291  #Jet_leptonPt = cms.string("userFloat('leptonPt')"),
292  Jet_leptonPtRel = cms.string("userFloat('leptonPtRelv0')"),
293  Jet_leptonPtRelInv = cms.string("userFloat('leptonPtRelInvv0')*jecFactor('Uncorrected')"),
294  Jet_leptonDeltaR = cms.string("userFloat('leptonDeltaR')"),
295  #Jet_leptonPdgId = cms.string("userInt('leptonPdgId')"),
296  Jet_neHEF = cms.string("neutralHadronEnergyFraction()"),
297  Jet_neEmEF = cms.string("neutralEmEnergyFraction()"),
298  Jet_chHEF = cms.string("chargedHadronEnergyFraction()"),
299  Jet_chEmEF = cms.string("chargedEmEnergyFraction()"),
300  isMu = cms.string("?abs(userInt('leptonPdgId'))==13?1:0"),
301  isEle = cms.string("?abs(userInt('leptonPdgId'))==11?1:0"),
302  isOther = cms.string("?userInt('leptonPdgId')==0?1:0"),
303  ),
304  inputTensorName = cms.string("ffwd_inp"),
305  outputTensorName = cms.string("ffwd_out/BiasAdd"),
306  outputNames = cms.vstring(["corr","res"]),
307  outputFormulas = cms.vstring(["at(0)*0.28492164611816406+1.0596693754196167","0.5*(at(2)-at(1))*0.28492164611816406"]),
308  nThreads = cms.uint32(1),
309  singleThreadPool = cms.string("no_threads"),
310 )
311 
312 
313 
314 ##### Soft Activity tables
315 saJetTable = cms.EDProducer("SimpleCandidateFlatTableProducer",
316  src = cms.InputTag("softActivityJets"),
317  cut = cms.string(""),
318  maxLen = cms.uint32(6),
319  name = cms.string("SoftActivityJet"),
320  doc = cms.string("jets clustered from charged candidates compatible with primary vertex (" + chsForSATkJets.cut.value()+")"),
321  singleton = cms.bool(False), # the number of entries is variable
322  extension = cms.bool(False), # this is the main table for the jets
323  variables = cms.PSet(P3Vars,
324  )
325 )
326 
327 saJetTable.variables.pt.precision=10
328 saJetTable.variables.eta.precision=8
329 saJetTable.variables.phi.precision=8
330 
331 saTable = cms.EDProducer("GlobalVariablesTableProducer",
332  variables = cms.PSet(
333  SoftActivityJetHT = ExtVar( cms.InputTag("softActivityJets"), "candidatescalarsum", doc = "scalar sum of soft activity jet pt, pt>1" ),
334  SoftActivityJetHT10 = ExtVar( cms.InputTag("softActivityJets10"), "candidatescalarsum", doc = "scalar sum of soft activity jet pt , pt >10" ),
335  SoftActivityJetHT5 = ExtVar( cms.InputTag("softActivityJets5"), "candidatescalarsum", doc = "scalar sum of soft activity jet pt, pt>5" ),
336  SoftActivityJetHT2 = ExtVar( cms.InputTag("softActivityJets2"), "candidatescalarsum", doc = "scalar sum of soft activity jet pt, pt >2" ),
337  SoftActivityJetNjets10 = ExtVar( cms.InputTag("softActivityJets10"), "candidatesize", doc = "number of soft activity jet pt, pt >2" ),
338  SoftActivityJetNjets5 = ExtVar( cms.InputTag("softActivityJets5"), "candidatesize", doc = "number of soft activity jet pt, pt >5" ),
339  SoftActivityJetNjets2 = ExtVar( cms.InputTag("softActivityJets2"), "candidatesize", doc = "number of soft activity jet pt, pt >10" ),
340 
341  )
342 )
343 
344 
345 
346 ## BOOSTED STUFF #################
347 fatJetTable = cms.EDProducer("SimpleCandidateFlatTableProducer",
348  src = cms.InputTag("finalJetsAK8"),
349  cut = cms.string(" pt > 170"), #probably already applied in miniaod
350  name = cms.string("FatJet"),
351  doc = cms.string("slimmedJetsAK8, i.e. ak8 fat jets for boosted analysis"),
352  singleton = cms.bool(False), # the number of entries is variable
353  extension = cms.bool(False), # this is the main table for the jets
354  variables = cms.PSet(P4Vars,
355  jetId = Var("userInt('tightId')*2+4*userInt('tightIdLepVeto')",int,doc="Jet ID flags bit1 is loose (always false in 2017 since it does not exist), bit2 is tight, bit3 is tightLepVeto"),
356  area = Var("jetArea()", float, doc="jet catchment area, for JECs",precision=10),
357  rawFactor = Var("1.-jecFactor('Uncorrected')",float,doc="1 - Factor to get back to raw pT",precision=6),
358  tau1 = Var("userFloat('NjettinessAK8Puppi:tau1')",float, doc="Nsubjettiness (1 axis)",precision=10),
359  tau2 = Var("userFloat('NjettinessAK8Puppi:tau2')",float, doc="Nsubjettiness (2 axis)",precision=10),
360  tau3 = Var("userFloat('NjettinessAK8Puppi:tau3')",float, doc="Nsubjettiness (3 axis)",precision=10),
361  tau4 = Var("userFloat('NjettinessAK8Puppi:tau4')",float, doc="Nsubjettiness (4 axis)",precision=10),
362  n2b1 = Var("userFloat('ak8PFJetsPuppiSoftDropValueMap:nb1AK8PuppiSoftDropN2')", float, doc="N2 with beta=1", precision=10),
363  n3b1 = Var("userFloat('ak8PFJetsPuppiSoftDropValueMap:nb1AK8PuppiSoftDropN3')", float, doc="N3 with beta=1", precision=10),
364  msoftdrop = Var("groomedMass('SoftDropPuppi')",float, doc="Corrected soft drop mass with PUPPI",precision=10),
365  btagCMVA = Var("bDiscriminator('pfCombinedMVAV2BJetTags')",float,doc="CMVA V2 btag discriminator",precision=10),
366  btagDeepB = Var("bDiscriminator('pfDeepCSVJetTags:probb')+bDiscriminator('pfDeepCSVJetTags:probbb')",float,doc="DeepCSV b+bb tag discriminator",precision=10),
367  btagCSVV2 = Var("bDiscriminator('pfCombinedInclusiveSecondaryVertexV2BJetTags')",float,doc=" pfCombinedInclusiveSecondaryVertexV2 b-tag discriminator (aka CSVV2)",precision=10),
368  btagHbb = Var("bDiscriminator('pfBoostedDoubleSecondaryVertexAK8BJetTags')",float,doc="Higgs to BB tagger discriminator",precision=10),
369  btagDDBvL_noMD = Var("bDiscriminator('pfDeepDoubleBvLJetTags:probHbb')",float,doc="DeepDoubleX discriminator (no mass-decorrelation) for H(Z)->bb vs QCD",precision=10),
370  btagDDCvL_noMD = Var("bDiscriminator('pfDeepDoubleCvLJetTags:probHcc')",float,doc="DeepDoubleX discriminator (no mass-decorrelation) for H(Z)->cc vs QCD",precision=10),
371  btagDDCvB_noMD = Var("bDiscriminator('pfDeepDoubleCvBJetTags:probHcc')",float,doc="DeepDoubleX discriminator (no mass-decorrelation) for H(Z)->cc vs H(Z)->bb",precision=10),
372  btagDDBvL = Var("bDiscriminator('pfMassIndependentDeepDoubleBvLJetTags:probHbb')",float,doc="DeepDoubleX (mass-decorrelated) discriminator for H(Z)->bb vs QCD",precision=10),
373  btagDDCvL = Var("bDiscriminator('pfMassIndependentDeepDoubleCvLJetTags:probHcc')",float,doc="DeepDoubleX (mass-decorrelated) discriminator for H(Z)->cc vs QCD",precision=10),
374  btagDDCvB = Var("bDiscriminator('pfMassIndependentDeepDoubleCvBJetTags:probHcc')",float,doc="DeepDoubleX (mass-decorrelated) discriminator for H(Z)->cc vs H(Z)->bb",precision=10),
375  deepTag_TvsQCD = Var("bDiscriminator('pfDeepBoostedDiscriminatorsJetTags:TvsQCD')",float,doc="DeepBoostedJet tagger top vs QCD discriminator",precision=10),
376  deepTag_WvsQCD = Var("bDiscriminator('pfDeepBoostedDiscriminatorsJetTags:WvsQCD')",float,doc="DeepBoostedJet tagger W vs QCD discriminator",precision=10),
377  deepTag_ZvsQCD = Var("bDiscriminator('pfDeepBoostedDiscriminatorsJetTags:ZvsQCD')",float,doc="DeepBoostedJet tagger Z vs QCD discriminator",precision=10),
378  deepTag_H = Var("bDiscriminator('pfDeepBoostedJetTags:probHbb')+bDiscriminator('pfDeepBoostedJetTags:probHcc')+bDiscriminator('pfDeepBoostedJetTags:probHqqqq')",float,doc="DeepBoostedJet tagger H(bb,cc,4q) sum",precision=10),
379  deepTag_QCD = Var("bDiscriminator('pfDeepBoostedJetTags:probQCDbb')+bDiscriminator('pfDeepBoostedJetTags:probQCDcc')+bDiscriminator('pfDeepBoostedJetTags:probQCDb')+bDiscriminator('pfDeepBoostedJetTags:probQCDc')+bDiscriminator('pfDeepBoostedJetTags:probQCDothers')",float,doc="DeepBoostedJet tagger QCD(bb,cc,b,c,others) sum",precision=10),
380  deepTag_QCDothers = Var("bDiscriminator('pfDeepBoostedJetTags:probQCDothers')",float,doc="DeepBoostedJet tagger QCDothers value",precision=10),
381  deepTagMD_TvsQCD = Var("bDiscriminator('pfMassDecorrelatedDeepBoostedDiscriminatorsJetTags:TvsQCD')",float,doc="Mass-decorrelated DeepBoostedJet tagger top vs QCD discriminator",precision=10),
382  deepTagMD_WvsQCD = Var("bDiscriminator('pfMassDecorrelatedDeepBoostedDiscriminatorsJetTags:WvsQCD')",float,doc="Mass-decorrelated DeepBoostedJet tagger W vs QCD discriminator",precision=10),
383  deepTagMD_ZvsQCD = Var("bDiscriminator('pfMassDecorrelatedDeepBoostedDiscriminatorsJetTags:ZvsQCD')",float,doc="Mass-decorrelated DeepBoostedJet tagger Z vs QCD discriminator",precision=10),
384  deepTagMD_ZHbbvsQCD = Var("bDiscriminator('pfMassDecorrelatedDeepBoostedDiscriminatorsJetTags:ZHbbvsQCD')",float,doc="Mass-decorrelated DeepBoostedJet tagger Z/H->bb vs QCD discriminator",precision=10),
385  deepTagMD_ZbbvsQCD = Var("bDiscriminator('pfMassDecorrelatedDeepBoostedDiscriminatorsJetTags:ZbbvsQCD')",float,doc="Mass-decorrelated DeepBoostedJet tagger Z->bb vs QCD discriminator",precision=10),
386  deepTagMD_HbbvsQCD = Var("bDiscriminator('pfMassDecorrelatedDeepBoostedDiscriminatorsJetTags:HbbvsQCD')",float,doc="Mass-decorrelated DeepBoostedJet tagger H->bb vs QCD discriminator",precision=10),
387  deepTagMD_ZHccvsQCD = Var("bDiscriminator('pfMassDecorrelatedDeepBoostedDiscriminatorsJetTags:ZHccvsQCD')",float,doc="Mass-decorrelated DeepBoostedJet tagger Z/H->cc vs QCD discriminator",precision=10),
388  deepTagMD_H4qvsQCD = Var("bDiscriminator('pfMassDecorrelatedDeepBoostedDiscriminatorsJetTags:H4qvsQCD')",float,doc="Mass-decorrelated DeepBoostedJet tagger H->4q vs QCD discriminator",precision=10),
389  deepTagMD_bbvsLight = Var("bDiscriminator('pfMassDecorrelatedDeepBoostedDiscriminatorsJetTags:bbvsLight')",float,doc="Mass-decorrelated DeepBoostedJet tagger Z/H/gluon->bb vs light flavour discriminator",precision=10),
390  deepTagMD_ccvsLight = Var("bDiscriminator('pfMassDecorrelatedDeepBoostedDiscriminatorsJetTags:ccvsLight')",float,doc="Mass-decorrelated DeepBoostedJet tagger Z/H/gluon->cc vs light flavour discriminator",precision=10),
391  subJetIdx1 = Var("?nSubjetCollections()>0 && subjets('SoftDropPuppi').size()>0?subjets('SoftDropPuppi')[0].key():-1", int,
392  doc="index of first subjet"),
393  subJetIdx2 = Var("?nSubjetCollections()>0 && subjets('SoftDropPuppi').size()>1?subjets('SoftDropPuppi')[1].key():-1", int,
394  doc="index of second subjet"),
395 
396 # btagDeepC = Var("bDiscriminator('pfDeepCSVJetTags:probc')",float,doc="CMVA V2 btag discriminator",precision=10),
397 #puIdDisc = Var("userFloat('pileupJetId:fullDiscriminant')",float,doc="Pilup ID discriminant",precision=10),
398 # nConstituents = Var("numberOfDaughters()",int,doc="Number of particles in the jet"),
399 # rawFactor = Var("1.-jecFactor('Uncorrected')",float,doc="1 - Factor to get back to raw pT",precision=6),
400  )
401 )
402 ### Era dependent customization
403 run2_miniAOD_80XLegacy.toModify( bjetNN,outputFormulas = cms.vstring(["at(0)*0.31628304719924927+1.0454729795455933","0.5*(at(2)-at(1))*0.31628304719924927"]))
404 run2_miniAOD_80XLegacy.toModify( fatJetTable.variables, msoftdrop_chs = Var("userFloat('ak8PFJetsCHSSoftDropMass')",float, doc="Legacy uncorrected soft drop mass with CHS",precision=10))
405 run2_miniAOD_80XLegacy.toModify( fatJetTable.variables.tau1, expr = cms.string("userFloat(\'ak8PFJetsPuppiValueMap:NjettinessAK8PuppiTau1\')"),)
406 run2_miniAOD_80XLegacy.toModify( fatJetTable.variables.tau2, expr = cms.string("userFloat(\'ak8PFJetsPuppiValueMap:NjettinessAK8PuppiTau2\')"),)
407 run2_miniAOD_80XLegacy.toModify( fatJetTable.variables.tau3, expr = cms.string("userFloat(\'ak8PFJetsPuppiValueMap:NjettinessAK8PuppiTau3\')"),)
408 run2_miniAOD_80XLegacy.toModify( fatJetTable.variables, tau4 = None)
409 run2_miniAOD_80XLegacy.toModify( fatJetTable.variables, n2b1 = None)
410 run2_miniAOD_80XLegacy.toModify( fatJetTable.variables, n3b1 = None)
411 run2_jme_2016.toModify( fatJetTable.variables, jetId = Var("userInt('tightId')*2+userInt('looseId')",int,doc="Jet ID flags bit1 is loose, bit2 is tight"))
412 
413 
414 
415 
416 subJetTable = cms.EDProducer("SimpleCandidateFlatTableProducer",
417  src = cms.InputTag("slimmedJetsAK8PFPuppiSoftDropPacked","SubJets"),
418  cut = cms.string(""), #probably already applied in miniaod
419  name = cms.string("SubJet"),
420  doc = cms.string("slimmedJetsAK8, i.e. ak8 fat jets for boosted analysis"),
421  singleton = cms.bool(False), # the number of entries is variable
422  extension = cms.bool(False), # this is the main table for the jets
423  variables = cms.PSet(P4Vars,
424  btagCMVA = Var("bDiscriminator('pfCombinedMVAV2BJetTags')",float,doc="CMVA V2 btag discriminator",precision=10),
425  btagDeepB = Var("bDiscriminator('pfDeepCSVJetTags:probb')+bDiscriminator('pfDeepCSVJetTags:probbb')",float,doc="DeepCSV b+bb tag discriminator",precision=10),
426  btagCSVV2 = Var("bDiscriminator('pfCombinedInclusiveSecondaryVertexV2BJetTags')",float,doc=" pfCombinedInclusiveSecondaryVertexV2 b-tag discriminator (aka CSVV2)",precision=10),
427  rawFactor = Var("1.-jecFactor('Uncorrected')",float,doc="1 - Factor to get back to raw pT",precision=6),
428  tau1 = Var("userFloat('NjettinessAK8Subjets:tau1')",float, doc="Nsubjettiness (1 axis)",precision=10),
429  tau2 = Var("userFloat('NjettinessAK8Subjets:tau2')",float, doc="Nsubjettiness (2 axis)",precision=10),
430  tau3 = Var("userFloat('NjettinessAK8Subjets:tau3')",float, doc="Nsubjettiness (3 axis)",precision=10),
431  tau4 = Var("userFloat('NjettinessAK8Subjets:tau4')",float, doc="Nsubjettiness (4 axis)",precision=10),
432  n2b1 = Var("userFloat('nb1AK8PuppiSoftDropSubjets:ecfN2')", float, doc="N2 with beta=1", precision=10),
433  n3b1 = Var("userFloat('nb1AK8PuppiSoftDropSubjets:ecfN3')", float, doc="N3 with beta=1", precision=10),
434  )
435 )
436 
437 #jets are not as precise as muons
438 fatJetTable.variables.pt.precision=10
439 subJetTable.variables.pt.precision=10
440 
441 run2_miniAOD_80XLegacy.toModify( subJetTable.variables, tau1 = None)
442 run2_miniAOD_80XLegacy.toModify( subJetTable.variables, tau2 = None)
443 run2_miniAOD_80XLegacy.toModify( subJetTable.variables, tau3 = None)
444 run2_miniAOD_80XLegacy.toModify( subJetTable.variables, tau4 = None)
445 run2_miniAOD_80XLegacy.toModify( subJetTable.variables, n2b1 = None)
446 run2_miniAOD_80XLegacy.toModify( subJetTable.variables, n3b1 = None)
447 run2_miniAOD_80XLegacy.toModify( subJetTable.variables, btagCMVA = None, btagDeepB = None)
448 
449 
450 corrT1METJetTable = cms.EDProducer("SimpleCandidateFlatTableProducer",
451  src = cms.InputTag("corrT1METJets"),
452  cut = cms.string(""),
453  name = cms.string("CorrT1METJet"),
454  doc = cms.string("Additional low-pt jets for Type-1 MET re-correction"),
455  singleton = cms.bool(False), # the number of entries is variable
456  extension = cms.bool(False), # this is the main table for the jets
457  variables = cms.PSet(
458  rawPt = Var("pt()*jecFactor('Uncorrected')",float,precision=10),
459  eta = Var("eta", float,precision=12),
460  phi = Var("phi", float, precision=12),
461  area = Var("jetArea()", float, doc="jet catchment area, for JECs",precision=10),
462  )
463 )
464 
465 
466 
467 ## MC STUFF ######################
468 jetMCTable = cms.EDProducer("SimpleCandidateFlatTableProducer",
469  src = cms.InputTag("linkedObjects","jets"),
470  cut = cms.string(""), #we should not filter on cross linked collections
471  name = cms.string("Jet"),
472  singleton = cms.bool(False), # the number of entries is variable
473  extension = cms.bool(True), # this is an extension table for the jets
474  variables = cms.PSet(
475  partonFlavour = Var("partonFlavour()", int, doc="flavour from parton matching"),
476  hadronFlavour = Var("hadronFlavour()", int, doc="flavour from hadron ghost clustering"),
477  genJetIdx = Var("?genJetFwdRef().backRef().isNonnull()?genJetFwdRef().backRef().key():-1", int, doc="index of matched gen jet"),
478  )
479 )
480 genJetTable = cms.EDProducer("SimpleCandidateFlatTableProducer",
481  src = cms.InputTag("slimmedGenJets"),
482  cut = cms.string("pt > 10"),
483  name = cms.string("GenJet"),
484  doc = cms.string("slimmedGenJets, i.e. ak4 Jets made with visible genparticles"),
485  singleton = cms.bool(False), # the number of entries is variable
486  extension = cms.bool(False), # this is the main table for the genjets
487  variables = cms.PSet(P4Vars,
488  #anything else?
489  )
490 )
491 patJetPartons = cms.EDProducer('HadronAndPartonSelector',
492  src = cms.InputTag("generator"),
493  particles = cms.InputTag("prunedGenParticles"),
494  partonMode = cms.string("Auto"),
495  fullChainPhysPartons = cms.bool(True)
496 )
497 genJetFlavourAssociation = cms.EDProducer("JetFlavourClustering",
498  jets = genJetTable.src,
499  bHadrons = cms.InputTag("patJetPartons","bHadrons"),
500  cHadrons = cms.InputTag("patJetPartons","cHadrons"),
501  partons = cms.InputTag("patJetPartons","physicsPartons"),
502  leptons = cms.InputTag("patJetPartons","leptons"),
503  jetAlgorithm = cms.string("AntiKt"),
504  rParam = cms.double(0.4),
505  ghostRescaling = cms.double(1e-18),
506  hadronFlavourHasPriority = cms.bool(False)
507 )
508 genJetFlavourTable = cms.EDProducer("GenJetFlavourTableProducer",
509  name = genJetTable.name,
510  src = genJetTable.src,
511  cut = genJetTable.cut,
512  deltaR = cms.double(0.1),
513  jetFlavourInfos = cms.InputTag("slimmedGenJetsFlavourInfos"),
514 )
515 
516 genJetAK8Table = cms.EDProducer("SimpleCandidateFlatTableProducer",
517  src = cms.InputTag("slimmedGenJetsAK8"),
518  cut = cms.string("pt > 100."),
519  name = cms.string("GenJetAK8"),
520  doc = cms.string("slimmedGenJetsAK8, i.e. ak8 Jets made with visible genparticles"),
521  singleton = cms.bool(False), # the number of entries is variable
522  extension = cms.bool(False), # this is the main table for the genjets
523  variables = cms.PSet(P4Vars,
524  #anything else?
525  )
526 )
527 genJetAK8FlavourAssociation = cms.EDProducer("JetFlavourClustering",
528  jets = genJetAK8Table.src,
529  bHadrons = cms.InputTag("patJetPartons","bHadrons"),
530  cHadrons = cms.InputTag("patJetPartons","cHadrons"),
531  partons = cms.InputTag("patJetPartons","physicsPartons"),
532  leptons = cms.InputTag("patJetPartons","leptons"),
533  jetAlgorithm = cms.string("AntiKt"),
534  rParam = cms.double(0.8),
535  ghostRescaling = cms.double(1e-18),
536  hadronFlavourHasPriority = cms.bool(False)
537 )
538 genJetAK8FlavourTable = cms.EDProducer("GenJetFlavourTableProducer",
539  name = genJetAK8Table.name,
540  src = genJetAK8Table.src,
541  cut = genJetAK8Table.cut,
542  deltaR = cms.double(0.1),
543  jetFlavourInfos = cms.InputTag("genJetAK8FlavourAssociation"),
544 )
545 genSubJetAK8Table = cms.EDProducer("SimpleCandidateFlatTableProducer",
546  src = cms.InputTag("slimmedGenJetsAK8SoftDropSubJets"),
547  cut = cms.string(""), ## These don't get a pt cut, but in miniAOD only subjets from fat jets with pt > 100 are kept
548  name = cms.string("SubGenJetAK8"),
549  doc = cms.string("slimmedGenJetsAK8SoftDropSubJets, i.e. subjets of ak8 Jets made with visible genparticles"),
550  singleton = cms.bool(False), # the number of entries is variable
551  extension = cms.bool(False), # this is the main table for the genjets
552  variables = cms.PSet(P4Vars,
553  #anything else?
554  )
555 )
556 ### Era dependent customization
557 run2_miniAOD_80XLegacy.toModify( genJetFlavourTable, jetFlavourInfos = cms.InputTag("genJetFlavourAssociation"),)
558 
559 from RecoJets.JetProducers.QGTagger_cfi import QGTagger
560 qgtagger=QGTagger.clone(srcJets="updatedJets",srcVertexCollection="offlineSlimmedPrimaryVertices")
561 
562 #before cross linking
563 jetSequence = cms.Sequence(jetCorrFactorsNano+updatedJets+tightJetId+tightJetIdLepVeto+bJetVars+jercVars+qgtagger+updatedJetsWithUserData+jetCorrFactorsAK8+updatedJetsAK8+tightJetIdAK8+tightJetIdLepVetoAK8+updatedJetsAK8WithUserData+chsForSATkJets+softActivityJets+softActivityJets2+softActivityJets5+softActivityJets10+finalJets+finalJetsAK8)
564 
565 _jetSequence_2016 = jetSequence.copy()
566 _jetSequence_2016.insert(_jetSequence_2016.index(tightJetId), looseJetId)
567 _jetSequence_2016.insert(_jetSequence_2016.index(tightJetIdAK8), looseJetIdAK8)
568 run2_jme_2016.toReplaceWith(jetSequence, _jetSequence_2016)
569 
570 #after cross linkining
571 jetTables = cms.Sequence(bjetMVA+bjetNN+jetTable+fatJetTable+subJetTable+saJetTable+saTable)
572 
573 #MC only producers and tables
574 jetMC = cms.Sequence(jetMCTable+genJetTable+patJetPartons+genJetFlavourTable+genJetAK8Table+genJetAK8FlavourAssociation+genJetAK8FlavourTable+genSubJetAK8Table)
575 _jetMC_pre94X = jetMC.copy()
576 _jetMC_pre94X.insert(_jetMC_pre94X.index(genJetFlavourTable),genJetFlavourAssociation)
577 _jetMC_pre94X.remove(genSubJetAK8Table)
578 run2_miniAOD_80XLegacy.toReplaceWith(jetMC, _jetMC_pre94X)
579 
580 
def ExtVar(tag, valtype, compression=None, doc=None, mcOnly=False, precision=-1)
Definition: common_cff.py:31
def Var(expr, valtype, compression=None, doc=None, mcOnly=False, precision=-1)
Definition: common_cff.py:20