Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

BTV fix for Puppi integration #4

Merged
merged 1 commit into from
Jan 24, 2020
Merged

BTV fix for Puppi integration #4

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
22 changes: 21 additions & 1 deletion PhysicsTools/JetMCAlgos/plugins/JetFlavourClustering.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -95,6 +95,7 @@
#include "DataFormats/HepMCCandidate/interface/GenParticleFwd.h"
#include "DataFormats/Math/interface/deltaR.h"
#include "PhysicsTools/JetMCUtils/interface/CandMCTag.h"
#include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"

#include "fastjet/JetDefinition.hh"
#include "fastjet/ClusterSequence.hh"
Expand Down Expand Up @@ -196,6 +197,8 @@ class JetFlavourClustering : public edm::stream::EDProducer<> {
const double relPtTolerance_;
const bool hadronFlavourHasPriority_;
const bool useSubjets_;
bool usePuppi_;

const bool useLeptons_;

ClusterSequencePtr fjClusterSeq_;
Expand All @@ -218,6 +221,7 @@ JetFlavourClustering::JetFlavourClustering(const edm::ParameterSet& iConfig)
partonsToken_(consumes<reco::GenParticleRefVector>(iConfig.getParameter<edm::InputTag>("partons"))),
jetAlgorithm_(iConfig.getParameter<std::string>("jetAlgorithm")),
rParam_(iConfig.getParameter<double>("rParam")),

jetPtMin_(
0.), // hardcoded to 0. since we simply want to recluster all input jets which already had some PtMin applied
ghostRescaling_(iConfig.exists("ghostRescaling") ? iConfig.getParameter<double>("ghostRescaling") : 1e-18),
Expand All @@ -232,6 +236,13 @@ JetFlavourClustering::JetFlavourClustering(const edm::ParameterSet& iConfig)
{
// register your products
produces<reco::JetFlavourInfoMatchingCollection>();
std::string label = iConfig.getParameter<edm::InputTag>("jets").label();
usePuppi_ = false;
if(label.find("Puppi") != std::string::npos){
usePuppi_ = true;
// This check will not fire on updatedPatJetsSlimmedDeepFlavour, updatedPatJetsSlimmedAK8DeepTags. Is that what we want?
}

if (useSubjets_)
produces<reco::JetFlavourInfoMatchingCollection>("SubJets");

Expand Down Expand Up @@ -289,6 +300,7 @@ void JetFlavourClustering::produce(edm::Event& iEvent, const edm::EventSetup& iS
if (useLeptons_)
iEvent.getByToken(leptonsToken_, leptons);


auto jetFlavourInfos = std::make_unique<reco::JetFlavourInfoMatchingCollection>(reco::JetRefBaseProd(jets));
std::unique_ptr<reco::JetFlavourInfoMatchingCollection> subjetFlavourInfos;
if (useSubjets_)
Expand All @@ -315,7 +327,15 @@ void JetFlavourClustering::produce(edm::Event& iEvent, const edm::EventSetup& iS
edm::LogWarning("NullTransverseMomentum") << "dropping input candidate with pt=0";
continue;
}
fjInputs.push_back(fastjet::PseudoJet(constit->px(), constit->py(), constit->pz(), constit->energy()));
if(usePuppi_){
const auto pf_constit = dynamic_cast<const reco::PFCandidate*>(&*constit);
double w = pf_constit->puppiWeight();
double E_w = std::sqrt(pf_constit->p() * w * pf_constit->p() * w + pf_constit->mass() * pf_constit->mass());
fjInputs.push_back(fastjet::PseudoJet(pf_constit->px() * w, pf_constit->py() * w, pf_constit->pz()* w, E_w));
}
else{
fjInputs.push_back(fastjet::PseudoJet(constit->px(), constit->py(), constit->pz(), constit->energy()));
}
}
}
// insert "ghost" b hadrons in the vector of constituents
Expand Down
64 changes: 53 additions & 11 deletions RecoBTag/SecondaryVertex/plugins/TemplatedSecondaryVertexProducer.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,11 @@
#include "FWCore/MessageLogger/interface/MessageLogger.h"
#include "FWCore/Utilities/interface/Exception.h"

#include "FWCore/ParameterSet/interface/Registry.h"
#include "FWCore/Common/interface/Provenance.h"
#include "DataFormats/Provenance/interface/ProductProvenance.h"


#include "DataFormats/PatCandidates/interface/Jet.h"
#include "DataFormats/TrackReco/interface/Track.h"
#include "DataFormats/TrackReco/interface/TrackFwd.h"
Expand Down Expand Up @@ -159,6 +164,7 @@ class TemplatedSecondaryVertexProducer : public edm::stream::EDProducer<> {
double ghostRescaling;
double relPtTolerance;
bool useFatJets;
bool usePuppi;
bool useGroomedFatJets;
edm::EDGetTokenT<edm::View<reco::Jet> > token_fatJets;
edm::EDGetTokenT<edm::View<reco::Jet> > token_groomedFatJets;
Expand Down Expand Up @@ -261,6 +267,7 @@ TemplatedSecondaryVertexProducer<IPTI, VTX>::TemplatedSecondaryVertexProducer(co
}
useSVClustering = (params.existsAs<bool>("useSVClustering") ? params.getParameter<bool>("useSVClustering") : false);
useSVMomentum = (params.existsAs<bool>("useSVMomentum") ? params.getParameter<bool>("useSVMomentum") : false);
usePuppi = false;
useFatJets = (useExternalSV && params.exists("fatJets"));
useGroomedFatJets = (useExternalSV && params.exists("groomedFatJets"));
if (useSVClustering) {
Expand All @@ -286,10 +293,16 @@ TemplatedSecondaryVertexProducer<IPTI, VTX>::TemplatedSecondaryVertexProducer(co
throw cms::Exception("InvalidJetAlgorithm") << "Jet clustering algorithm is invalid: " << jetAlgorithm
<< ", use CambridgeAachen | Kt | AntiKt" << std::endl;
}
if (useFatJets)
if (useFatJets){
token_fatJets = consumes<edm::View<reco::Jet> >(params.getParameter<edm::InputTag>("fatJets"));
if (useGroomedFatJets)
std::string label = params.getParameter<edm::InputTag>("fatJets").label();
if( label.find("Puppi") != std::string::npos ){
usePuppi = true;
}
}
if (useGroomedFatJets){
token_groomedFatJets = consumes<edm::View<reco::Jet> >(params.getParameter<edm::InputTag>("groomedFatJets"));
}
if (useFatJets && !useSVClustering)
rParam = params.getParameter<double>("rParam"); // will be used later as a dR cut

Expand All @@ -310,7 +323,21 @@ void TemplatedSecondaryVertexProducer<IPTI, VTX>::produce(edm::Event &event, con

edm::Handle<std::vector<IPTI> > trackIPTagInfos;
event.getByToken(token_trackIPTagInfo, trackIPTagInfos);

if(!useFatJets){
if constexpr (std::is_same_v<IPTI, CandIPTagInfo> ){
const edm::Provenance *prov = trackIPTagInfos.provenance();
const edm::ParameterSet& psetFromProvenance = edm::parameterSet(*prov,event.processHistory());
std::string label;
label = psetFromProvenance.getParameter<edm::InputTag>("jets").label();
if( label.find("Puppi") != std::string::npos ){
usePuppi = true;
}
}
else{
usePuppi = false;
}
}

// External Sec Vertex collection (e.g. for IVF usage)
edm::Handle<edm::View<VTX> > extSecVertex;
if (useExternalSV)
Expand All @@ -320,7 +347,6 @@ void TemplatedSecondaryVertexProducer<IPTI, VTX>::produce(edm::Event &event, con
edm::Handle<edm::View<reco::Jet> > groomedFatJetsHandle;
if (useFatJets) {
event.getByToken(token_fatJets, fatJetsHandle);

if (useGroomedFatJets) {
event.getByToken(token_groomedFatJets, groomedFatJetsHandle);

Expand Down Expand Up @@ -363,7 +389,7 @@ void TemplatedSecondaryVertexProducer<IPTI, VTX>::produce(edm::Event &event, con
default:
/* nothing */;
}

// ------------------------------------ SV clustering START --------------------------------------------
std::vector<std::vector<int> > clusteredSVs(trackIPTagInfos->size(), std::vector<int>());
if (useExternalSV && useSVClustering && !trackIPTagInfos->empty()) {
Expand All @@ -376,11 +402,19 @@ void TemplatedSecondaryVertexProducer<IPTI, VTX>::produce(edm::Event &event, con
std::vector<edm::Ptr<reco::Candidate> >::const_iterator m;
for (m = constituents.begin(); m != constituents.end(); ++m) {
reco::CandidatePtr constit = *m;
if (constit->pt() == 0) {
edm::LogWarning("NullTransverseMomentum") << "dropping input candidate with pt=0";
continue;
}
fjInputs.push_back(fastjet::PseudoJet(constit->px(), constit->py(), constit->pz(), constit->energy()));
if (constit->pt() == 0) {
edm::LogWarning("NullTransverseMomentum") << "dropping input candidate with pt=0";
continue;
}
if(usePuppi){
const reco::PFCandidate* pf_constit = dynamic_cast<const reco::PFCandidate*>(&*constit);
auto w = pf_constit->puppiWeight();
double E_w = std::sqrt(pf_constit->p() * w * pf_constit->p() * w + pf_constit->mass() * pf_constit->mass());
fjInputs.push_back(fastjet::PseudoJet(pf_constit->px() * w, pf_constit->py() * w, pf_constit->pz()* w, E_w));
}
else{
fjInputs.push_back(fastjet::PseudoJet(constit->px(), constit->py(), constit->pz(), constit->energy()));
}
}
}
} else {
Expand All @@ -394,7 +428,15 @@ void TemplatedSecondaryVertexProducer<IPTI, VTX>::produce(edm::Event &event, con
edm::LogWarning("NullTransverseMomentum") << "dropping input candidate with pt=0";
continue;
}
fjInputs.push_back(fastjet::PseudoJet(constit->px(), constit->py(), constit->pz(), constit->energy()));
if(usePuppi){
const auto* pf_constit = dynamic_cast<const reco::PFCandidate*>(&*constit);
auto w = pf_constit->puppiWeight();
double E_w = std::sqrt(pf_constit->p() * w * pf_constit->p() * w + pf_constit->mass() * pf_constit->mass());
fjInputs.push_back(fastjet::PseudoJet(pf_constit->px() * w, pf_constit->py() * w, pf_constit->pz()* w, E_w));
}
else{
fjInputs.push_back(fastjet::PseudoJet(constit->px(), constit->py(), constit->pz(), constit->energy()));
}
}
}
}
Expand Down