diff --git a/plugins/JetConstituentTableProducer.cc b/plugins/JetConstituentTableProducer.cc index 1fe8f5d..ae323e4 100644 --- a/plugins/JetConstituentTableProducer.cc +++ b/plugins/JetConstituentTableProducer.cc @@ -48,19 +48,24 @@ class JetConstituentTableProducer : public edm::stream::EDProducer<> { //const std::string name_; const std::string name_; const std::string nameSV_; + const std::string nameMu_; const std::string idx_name_; const std::string idx_nameSV_; + const std::string idx_nameMu_; const bool readBtag_; + const bool addMuonTable_; const double jet_radius_; edm::EDGetTokenT> jet_token_; edm::EDGetTokenT vtx_token_; edm::EDGetTokenT cand_token_; edm::EDGetTokenT sv_token_; + edm::EDGetTokenT > muon_token_; edm::Handle vtxs_; edm::Handle cands_; edm::Handle svs_; + edm::Handle > muons_; edm::ESHandle track_builder_; edm::ESGetToken track_builder_token_; @@ -75,18 +80,23 @@ template< typename T> JetConstituentTableProducer::JetConstituentTableProducer(const edm::ParameterSet &iConfig) : name_(iConfig.getParameter("name")), nameSV_(iConfig.getParameter("nameSV")), + nameMu_(iConfig.getParameter("nameMu")), idx_name_(iConfig.getParameter("idx_name")), idx_nameSV_(iConfig.getParameter("idx_nameSV")), + idx_nameMu_(iConfig.getParameter("idx_nameMu")), readBtag_(iConfig.getParameter("readBtag")), + addMuonTable_(iConfig.getParameter("addMuonTable")), jet_radius_(iConfig.getParameter("jet_radius")), jet_token_(consumes>(iConfig.getParameter("jets"))), vtx_token_(consumes(iConfig.getParameter("vertices"))), cand_token_(consumes(iConfig.getParameter("candidates"))), sv_token_(consumes(iConfig.getParameter("secondary_vertices"))), + muon_token_(consumes >(iConfig.getParameter("muons"))), track_builder_token_( esConsumes(edm::ESInputTag("", "TransientTrackBuilder"))){ produces(name_); produces(nameSV_); + produces(nameMu_); produces>(); } @@ -98,24 +108,33 @@ void JetConstituentTableProducer::produce(edm::Event &iEvent, const edm::Even // elements in all these collections must have the same order! auto outCands = std::make_unique>(); auto outSVs = std::make_unique> (); - std::vector jetIdx_pf, jetIdx_sv, pfcandIdx, svIdx; + auto outMuons = std::make_unique>(); + std::vector jetIdx_pf, jetIdx_sv, jetIdx_mu, pfcandIdx, svIdx, muIdx; // PF Cands std::vector btagEtaRel, btagPtRatio, btagPParRatio, btagSip3dVal, btagSip3dSig, btagJetDistVal, cand_pt; // Secondary vertices std::vector sv_mass, sv_pt, sv_ntracks, sv_chi2, sv_normchi2, sv_dxy, sv_dxysig, sv_d3d, sv_d3dsig, sv_costhetasvpv; std::vector sv_ptrel, sv_phirel, sv_deltaR, sv_enratio; + // Muons + std::vector muon_pt, muon_eta, muon_phi, muon_ptrel; + std::vector muon_isGlobal, muon_chi2Tk, muon_chi2; + std::vector muon_nMuHit, muon_nMatched, muon_nTkHit, muon_nPixHit, muon_nOutHit; auto jets = iEvent.getHandle(jet_token_); iEvent.getByToken(vtx_token_, vtxs_); iEvent.getByToken(cand_token_, cands_); iEvent.getByToken(sv_token_, svs_); + iEvent.getByToken(muon_token_, muons_); if(readBtag_){ track_builder_ = iSetup.getHandle(track_builder_token_); } + + for (unsigned i_jet = 0; i_jet < jets->size(); ++i_jet) { const auto &jet = jets->at(i_jet); + //std::cout<<"jet "<::produce(edm::Event &iEvent, const edm::Even } } } // end jet loop + + + //Muons + uint idx_mu=0; + for (const auto &mu : *muons_) { + if(reco::deltaR2(mu, jet) < jet_radius_ * jet_radius_){ + outMuons->push_back(mu); + if (readBtag_ && addMuonTable_) { + jetIdx_mu.push_back(i_jet); + muIdx.push_back(idx_mu); + muon_pt.push_back(mu.pt()); + muon_eta.push_back(mu.eta()); + muon_phi.push_back(mu.phi()); + muon_ptrel.push_back(mu.pt()/jet.pt()); + + if(mu.isGlobalMuon()){ + muon_nMuHit.push_back(mu.outerTrack()->hitPattern().numberOfValidMuonHits()); + muon_nMatched.push_back(mu.numberOfMatchedStations()); + muon_nTkHit.push_back(mu.innerTrack()->hitPattern().numberOfValidHits()); + muon_nPixHit.push_back(mu.innerTrack()->hitPattern().numberOfValidPixelHits()); + muon_nOutHit.push_back(mu.innerTrack()->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_OUTER_HITS)); + muon_chi2.push_back(mu.globalTrack()->normalizedChi2()); + muon_chi2Tk.push_back(mu.innerTrack()->normalizedChi2()); + } + else{ + muon_nMuHit.push_back(-1); + muon_nMatched.push_back(-1); + muon_nTkHit.push_back(-1); + muon_nPixHit.push_back(-1); + muon_nOutHit.push_back(-1); + muon_chi2.push_back(-1); + muon_chi2Tk.push_back(-1); + + } + muon_isGlobal.push_back(mu.isGlobalMuon()); + + } + idx_mu++; + } + + } } auto candTable = std::make_unique(outCands->size(), name_, false); @@ -256,6 +316,28 @@ void JetConstituentTableProducer::produce(edm::Event &iEvent, const edm::Even iEvent.put(std::move(svTable), nameSV_); iEvent.put(std::move(outCands)); + + + // Muon table + auto muonTable = std::make_unique(outMuons->size(), nameMu_, false); + // We fill from here only stuff that cannot be created with the SimpleFlatTnameableProducer + if (readBtag_ && addMuonTable_) { + muonTable->addColumn("jetIdx", jetIdx_mu, "Index of the parent jet"); + muonTable->addColumn(idx_nameMu_, muIdx, "Index in the Muon list"); + muonTable->addColumn("pt", muon_pt, "pt", 20); + muonTable->addColumn("eta", muon_eta, "eta", 20); + muonTable->addColumn("phi", muon_phi, "phi", 20); + muonTable->addColumn("isGlobal", muon_isGlobal, "isGlobal", 20); + muonTable->addColumn("ptrel", muon_ptrel, "pt relative to parent jet", 20); + muonTable->addColumn("nMuHit", muon_nMuHit, "number of muon hits", 20); + muonTable->addColumn("nMatched", muon_nMatched, "number of matched stations", 20); + muonTable->addColumn("nTkHit", muon_nTkHit, "number of tracker hits", 20); + muonTable->addColumn("nPixHit", muon_nPixHit, "number of pixel hits", 20); + muonTable->addColumn("nOutHit", muon_nOutHit, "number of missing outer hits", 20); + muonTable->addColumn("chi2", muon_chi2, "chi2", 20); + muonTable->addColumn("chi2Tk", muon_chi2Tk, "chi2 inner track", 20); + } + iEvent.put(std::move(muonTable), nameMu_); } template< typename T> @@ -263,15 +345,19 @@ void JetConstituentTableProducer::fillDescriptions(edm::ConfigurationDescript edm::ParameterSetDescription desc; desc.add("name", "JetPFCands"); desc.add("nameSV", "JetSV"); + desc.add("nameMu", "JetMuons"); desc.add("idx_name", "candIdx"); desc.add("idx_nameSV", "svIdx"); + desc.add("idx_nameMu", "muIdx"); desc.add("jet_radius", true); desc.add("readBtag", true); + desc.add("addMuonTable", false); desc.add("jets", edm::InputTag("slimmedJetsAK8")); desc.add("vertices", edm::InputTag("offlineSlimmedPrimaryVertices")); desc.add("candidates", edm::InputTag("packedPFCandidates")); desc.add("secondary_vertices", edm::InputTag("slimmedSecondaryVertices")); - descriptions.addWithDefaultLabel(desc); + desc.add("muons", edm::InputTag("slimmedMuons")); + descriptions.addWithDefaultLabel(desc); } typedef JetConstituentTableProducer PatJetConstituentTableProducer; diff --git a/python/addBTV.py b/python/addBTV.py index d24e4f3..0eed8b3 100644 --- a/python/addBTV.py +++ b/python/addBTV.py @@ -299,6 +299,84 @@ def add_BTV(process, runOnMC=False, onlyAK4=False, onlyAK8=False, keepInputs=['D doc="DeepCSV light btag discriminator", precision=10), ) + + FatJetVars = cms.PSet( + uncorrpt=Var("?availableJECSets().size()>0?correctedJet('Uncorrected').pt():pt()", + float, + doc="Uncorrected pT", + precision=10), + residual=Var("?availableJECSets().size()>0 ? pt()/correctedJet('L3Absolute').pt() : 1. ", + float, + doc="residual", + precision=10), + jes=Var("?availableJECSets().size()>0 ? pt()/correctedJet('Uncorrected').pt() : 1.", + float, + doc="jes (jet substructure)", + precision=10), + CombIVF=Var("bDiscriminator('pfCombinedInclusiveSecondaryVertexV2BJetTags')", + float, + doc="combinedIVF", + precision=10), + DeepCSVBDisc=Var("? (bDiscriminator('pfDeepCSVJetTags:probb')+bDiscriminator('pfDeepCSVJetTags:probbb'))>=0 ? bDiscriminator('pfDeepCSVJetTags:probb')+bDiscriminator('pfDeepCSVJetTags:probbb') : -1", + float, + doc="DeepCSVBDisc", + precision=10), + btagDeepB_c=Var("bDiscriminator('pfDeepCSVJetTags:probc')", + float, + doc="DeepCSV c tag discriminator", + precision=10), + DeepCSVCvsLDisc=Var("? (bDiscriminator('pfDeepCSVJetTags:probc'))>=0 ? bDiscriminator('pfDeepCSVJetTags:probc')/(bDiscriminator('pfDeepCSVJetTags:probc')+bDiscriminator('pfDeepCSVJetTags:probudsg')) : -1", + float, + doc="DeepCSV C vs L discriminator", + precision=10), + DeepCSVCvsBDisc=Var("? (bDiscriminator('pfDeepCSVJetTags:probc'))>=0 ? bDiscriminator('pfDeepCSVJetTags:probc')/(bDiscriminator('pfDeepCSVJetTags:probc')+bDiscriminator('pfDeepCSVJetTags:probb')+bDiscriminator('pfDeepCSVJetTags:probbb')) : -1", + float, + doc="DeepCSV C vs B discriminator", + precision=10), + + DoubleSV=Var("bDiscriminator('doubleSVBJetTags')", + float, + doc="DoubleSV discriminator", + precision=10), + ) + + SubJetVars = cms.PSet( + uncorrpt=Var("?availableJECSets().size()>0?correctedJet('Uncorrected').pt():pt()", + float, + doc="Uncorrected pT", + precision=10), + residual=Var("?availableJECSets().size()>0 ? pt()/correctedJet('L3Absolute').pt() : 1. ", + float, + doc="residual", + precision=10), + jes=Var("?availableJECSets().size()>0 ? pt()/correctedJet('Uncorrected').pt() : 1.", + float, + doc="jes (jet substructure)", + precision=10), + CombIVF=Var("bDiscriminator('pfCombinedInclusiveSecondaryVertexV2BJetTags')", + float, + doc="combinedIVF", + precision=10), + DeepCSVBDisc=Var("? (bDiscriminator('pfDeepCSVJetTags:probb')+bDiscriminator('pfDeepCSVJetTags:probbb'))>=0 ? bDiscriminator('pfDeepCSVJetTags:probb')+bDiscriminator('pfDeepCSVJetTags:probbb') : -1", + float, + doc="DeepCSVBDisc", + precision=10), + btagDeepB_c=Var("bDiscriminator('pfDeepCSVJetTags:probc')", + float, + doc="DeepCSV c tag discriminator", + precision=10), + DeepCSVCvsLDisc=Var("? (bDiscriminator('pfDeepCSVJetTags:probc'))>=0 ? bDiscriminator('pfDeepCSVJetTags:probc')/(bDiscriminator('pfDeepCSVJetTags:probc')+bDiscriminator('pfDeepCSVJetTags:probudsg')) : -1", + float, + doc="DeepCSV C vs L discriminator", + precision=10), + DeepCSVCvsBDisc=Var("? (bDiscriminator('pfDeepCSVJetTags:probc'))>=0 ? bDiscriminator('pfDeepCSVJetTags:probc')/(bDiscriminator('pfDeepCSVJetTags:probc')+bDiscriminator('pfDeepCSVJetTags:probb')+bDiscriminator('pfDeepCSVJetTags:probbb')) : -1", + float, + doc="DeepCSV C vs B discriminator", + precision=10), + + ) + + # decouple these from CommonVars, not relevant for data HadronCountingVars = cms.PSet( @@ -348,6 +426,7 @@ def add_BTV(process, runOnMC=False, onlyAK4=False, onlyAK8=False, keepInputs=['D extension=cms.bool(True), # this is the extension table for FatJets variables=cms.PSet( CommonVars, + FatJetVars, #HadronCountingVars if runOnMC else cms.PSet(), # only necessary before 106x get_DDX_vars() if ('DDX' in keepInputs) else cms.PSet(), )) @@ -363,6 +442,7 @@ def add_BTV(process, runOnMC=False, onlyAK4=False, onlyAK8=False, keepInputs=['D extension=cms.bool(True), # this is the extension table for FatJets variables=cms.PSet( CommonVars, + SubJetVars, #HadronCountingVars if runOnMC else cms.PSet(), # only necessary before 106x )) diff --git a/python/addPFCands_cff.py b/python/addPFCands_cff.py index 1cb8602..594e6dc 100644 --- a/python/addPFCands_cff.py +++ b/python/addPFCands_cff.py @@ -62,7 +62,10 @@ def addPFCands(process, runOnMC=False, allPF = False, onlyAK4=False, onlyAK8=Fal idx_name = cms.string("pFCandsIdx"), nameSV = cms.string("FatJetSVs"), idx_nameSV = cms.string("sVIdx"), - ) + nameMu = cms.string("FatJetMuons"), + idx_nameMu = cms.string("MuIdx"), + addMuonTable = cms.bool(True), + ) process.customAK4ConstituentsTable = cms.EDProducer("PatJetConstituentTableProducer", candidates = candInput, jets = cms.InputTag("finalJetsPuppi"), # was finalJets before @@ -71,6 +74,7 @@ def addPFCands(process, runOnMC=False, allPF = False, onlyAK4=False, onlyAK8=Fal idx_name = cms.string("pFCandsIdx"), nameSV = cms.string("JetSVs"), idx_nameSV = cms.string("sVIdx"), + addMuonTable = cms.bool(False), ) if not allPF: process.customizedPFCandsTask.add(process.finalJetsConstituents) @@ -122,7 +126,11 @@ def addPFCands(process, runOnMC=False, allPF = False, onlyAK4=False, onlyAK8=Fal nameSV = cms.string("GenFatJetSVs"), idx_name = cms.string("pFCandsIdx"), idx_nameSV = cms.string("sVIdx"), - readBtag = cms.bool(False)) + nameMu = cms.string("GenFatJetMuons"), + idx_nameMu = cms.string("MuIdx"), + readBtag = cms.bool(False), + addMuonTable = cms.bool(False), + ) process.genAK4ConstituentsTable = cms.EDProducer("GenJetConstituentTableProducer", candidates = genCandInput, jets = cms.InputTag("genJetsAK4Constituents"), # Note: The name has "Constituents" in it, but these are the jets @@ -130,7 +138,11 @@ def addPFCands(process, runOnMC=False, allPF = False, onlyAK4=False, onlyAK8=Fal nameSV = cms.string("GenJetSVs"), idx_name = cms.string("pFCandsIdx"), idx_nameSV = cms.string("sVIdx"), - readBtag = cms.bool(False)) + nameMu = cms.string("GenJetMuons"), + idx_nameMu = cms.string("MuIdx"), + readBtag = cms.bool(False), + addMuonTable = cms.bool(False), + ) process.customizedPFCandsTask.add(process.genJetsAK4Constituents) #Note: For gen need to add jets to the process to keep pt cuts. process.customizedPFCandsTask.add(process.genJetsAK8Constituents) if not allPF: diff --git a/test/nano_mc_Run3_122X_NANO.py b/test/nano_mc_Run3_122X_NANO.py index ef85edc..424a8bb 100644 --- a/test/nano_mc_Run3_122X_NANO.py +++ b/test/nano_mc_Run3_122X_NANO.py @@ -23,13 +23,13 @@ process.load('Configuration.StandardSequences.FrontierConditions_GlobalTag_cff') process.maxEvents = cms.untracked.PSet( - input = cms.untracked.int32(-1), + input = cms.untracked.int32(1000), output = cms.optional.untracked.allowed(cms.int32,cms.PSet) ) # Input source process.source = cms.Source("PoolSource", - fileNames = cms.untracked.vstring('/store/mc/Run3Winter22MiniAOD/TTTo2L2Nu_CP5_13p6TeV_powheg-pythia8/MINIAODSIM/122X_mcRun3_2021_realistic_v9-v2/2550000/0d44f6e9-6961-4d60-b2c1-0e21c1249100.root'), + fileNames = cms.untracked.vstring('/store/mc/Run3Winter22MiniAOD/ZprimeToTT_M3000_W30_TuneCP2_13p6TeV-madgraph-pythiaMLM-pythia8/MINIAODSIM/122X_mcRun3_2021_realistic_v9-v2/60000/0aca854d-fcdb-4523-834d-201c36e4a6ca.root'), secondaryFileNames = cms.untracked.vstring() )