Skip to content

Commit

Permalink
Use original PackedCandidate and other changes.
Browse files Browse the repository at this point in the history 
 - Use original PackedCandidate instead of puppi scaled ones. (needs
cms-sw/cmssw#28035)
 - Use uncorr_pt to sort subjets.
 - Use calibrated hcalFrac. (11_0_X)
 - Add Top_bl matching def.
 - Return matched parton for QCD labels.
  • Loading branch information
@hqucms
hqucms committed Sep 28, 2019
1 parent 09d3802 commit bb65a79
Show file tree
Hide file tree
Showing 16 changed files with 394 additions and 117 deletions.
4 changes: 2 additions & 2 deletions FatJetHelpers/interface/FatJetMatching.h
View file
Expand Up @@ -48,7 +48,7 @@ class FatJetMatching {

enum FatJetLabel {
Invalid=0,
Top_all=10, Top_bcq, Top_bqq, Top_bc, Top_bq,
Top_all=10, Top_bcq, Top_bqq, Top_bc, Top_bq, Top_bele, Top_bmu, Top_btau,
W_all=20, W_cq, W_qq,
Z_all=30, Z_bb, Z_cc, Z_qq,
H_all=40, H_bb, H_cc, H_qqqq, H_tautau,
Expand All @@ -70,7 +70,7 @@ class FatJetMatching {
std::pair<FatJetLabel, const reco::GenParticle*> w_label(const pat::Jet *jet, const reco::GenParticle *parton, double distR);
std::pair<FatJetLabel, const reco::GenParticle*> z_label(const pat::Jet *jet, const reco::GenParticle *parton, double distR);
std::pair<FatJetLabel, const reco::GenParticle*> higgs_label(const pat::Jet *jet, const reco::GenParticle *parton, double distR);
std::pair<FatJetLabel, const reco::GenParticle*> qcd_label(const pat::Jet *jet);
std::pair<FatJetLabel, const reco::GenParticle*> qcd_label(const pat::Jet *jet, const reco::GenParticleCollection& genParticles, double distR);


private:
Expand Down
77 changes: 70 additions & 7 deletions FatJetHelpers/src/FatJetMatching.cc
View file
Expand Up @@ -217,7 +217,7 @@ std::pair<FatJetMatching::FatJetLabel, const reco::GenParticle*> FatJetMatching:
if (genParticles.size() != processed_.size())
throw std::logic_error("[FatJetMatching::flavor] Not all genParticles are processed!");

return qcd_label(jet);
return qcd_label(jet, genParticles, distR);

}

Expand Down Expand Up @@ -369,6 +369,47 @@ std::pair<FatJetMatching::FatJetLabel,const reco::GenParticle*> FatJetMatching::
return w_label(jet, w_from_top, distR);
}
}
} else {
// leptonic W
if (debug_){
using namespace std;
cout << "jet: " << jet->polarP4() << endl;
cout << "top: "; printGenParticleInfo(top, -1);
cout << "b: "; printGenParticleInfo(b_from_top, -1);
cout << "W: "; printGenParticleInfo(w_from_top, -1);
}

const reco::GenParticle* lep = nullptr;
for (unsigned i=0; i<w_from_top->numberOfDaughters(); ++i){
const auto *dau = dynamic_cast<const reco::GenParticle*>(w_from_top->daughter(i));
auto pdgid = std::abs(dau->pdgId());
if (pdgid == ParticleID::p_eminus || pdgid == ParticleID::p_muminus || pdgid == ParticleID::p_tauminus){
// use final version here!
lep = getFinal(dau); break;
}
}

if (!lep) throw std::logic_error("[FatJetMatching::top_label] Cannot find charged lepton from leptonic W decay!");

double dr_b = reco::deltaR(jet->p4(), b_from_top->p4());
double dr_l = reco::deltaR(jet->p4(), lep->p4());
if (debug_){
using namespace std;
cout << "deltaR(jet, b) : " << dr_b << endl;
cout << "deltaR(jet, l) : " << dr_l << endl;
cout << "pdgid(l) : " << lep->pdgId() << endl;
}

if (dr_b < distR && dr_l < distR){
auto pdgid = std::abs(lep->pdgId());
if (pdgid == ParticleID::p_eminus){
return std::make_pair(FatJetLabel::Top_bele, top);
} else if (pdgid == ParticleID::p_muminus){
return std::make_pair(FatJetLabel::Top_bmu, top);
} else if (pdgid == ParticleID::p_tauminus){
return std::make_pair(FatJetLabel::Top_btau, top);
}
}
}

return std::make_pair(FatJetLabel::Invalid, nullptr);
Expand Down Expand Up @@ -577,6 +618,7 @@ std::pair<FatJetMatching::FatJetLabel,const reco::GenParticle*> FatJetMatching::
}
if (taus.size() == 2){
// higgs -> tautau
// use first version or last version of the tau in dr?
double dr_tau1 = reco::deltaR(jet->p4(), taus.at(0)->p4());
double dr_tau2 = reco::deltaR(jet->p4(), taus.at(1)->p4());

Expand Down Expand Up @@ -610,20 +652,41 @@ std::pair<FatJetMatching::FatJetLabel,const reco::GenParticle*> FatJetMatching::

}

std::pair<FatJetMatching::FatJetLabel,const reco::GenParticle*> FatJetMatching::qcd_label(const pat::Jet* jet)
std::pair<FatJetMatching::FatJetLabel,const reco::GenParticle*> FatJetMatching::qcd_label(const pat::Jet* jet, const reco::GenParticleCollection& genParticles, double distR)
{

const reco::GenParticle *parton = nullptr;
double minDR = 999;
for (const auto &gp : genParticles){
if (gp.status() != 23) continue;
auto pdgid = std::abs(gp.pdgId());
if (!(pdgid<ParticleID::p_t || pdgid==ParticleID::p_g)) continue;
auto dr = reco::deltaR(gp, *jet);
if (dr<distR && dr<minDR){
minDR = dr;
parton = &gp;
}
}
if (debug_){
using namespace std;
if (parton){
cout << "parton"; printGenParticleInfo(parton, -1);
cout << "dr(jet, parton): " << minDR << endl;
}
}

auto n_bHadrons = jet->jetFlavourInfo().getbHadrons().size();
auto n_cHadrons = jet->jetFlavourInfo().getcHadrons().size();

if (n_bHadrons>=2) {
return std::make_pair(FatJetLabel::QCD_bb, nullptr);
return std::make_pair(FatJetLabel::QCD_bb, parton);
}else if (n_bHadrons==1){
return std::make_pair(FatJetLabel::QCD_b, nullptr);
return std::make_pair(FatJetLabel::QCD_b, parton);
}else if (n_cHadrons>=2){
return std::make_pair(FatJetLabel::QCD_cc, nullptr);
return std::make_pair(FatJetLabel::QCD_cc, parton);
}else if (n_cHadrons==1){
return std::make_pair(FatJetLabel::QCD_c, nullptr);
return std::make_pair(FatJetLabel::QCD_c, parton);
}

return std::make_pair(FatJetLabel::QCD_others, nullptr);
return std::make_pair(FatJetLabel::QCD_others, parton);
}
22 changes: 16 additions & 6 deletions NtupleCommons/interface/JetHelper.h
View file
Expand Up @@ -12,12 +12,14 @@
#include "DataFormats/PatCandidates/interface/Jet.h"
#include "DataFormats/PatCandidates/interface/PackedCandidate.h"

#include <map>

namespace deepntuples {

class JetHelper {
public:
JetHelper() {}
JetHelper(const pat::Jet *jet, bool has_puppi_weighted_daughters);
JetHelper(const pat::Jet *jet, const edm::Handle<reco::CandidateView> &pfcands);

virtual ~JetHelper() {}

Expand All @@ -27,12 +29,19 @@ class JetHelper {
// ------

// return jet constituents (PF candidates)
const std::vector<const pat::PackedCandidate*>& getJetConstituents() const { return daughters_; }
const std::vector<reco::CandidatePtr>& getJetConstituents() const { return daughters_; }
unsigned int numberOfDaughters() const { return daughters_.size(); }
bool hasPuppiWeightedDaughters() const { return has_puppi_weighted_daughters_; }
float getPuppiWeight(const reco::CandidatePtr &cand) const {
auto iter = puppi_wgt_cache_.find(cand.key());
if (iter == puppi_wgt_cache_.end()){
throw cms::Exception("[JetHelper::getPuppiWeight] Cannot get puppi wgt!");
}
return iter->second;
}

const pat::Jet& jet() const { return *jet_; }
const std::vector<const pat::Jet*>& getSubJets() const { return subjets_; }
const std::vector<const pat::Jet*>& getUncorrSubJets() const { return uncorr_subjets_; }

const reco::GenJet* genjetWithNu() const { return genjetWithNu_; }
const reco::GenJet* genjetWithNuSoftDrop() const { return genjetWithNuSoftDrop_; }
Expand All @@ -41,17 +50,18 @@ class JetHelper {


private:
void initializeConstituents();
void initializeConstituents(const edm::Handle<reco::CandidateView> &pfcands);


private:
// data members
const pat::Jet *jet_ = nullptr;
bool has_puppi_weighted_daughters_ = false;
const reco::GenJet *genjetWithNu_ = nullptr;
const reco::GenJet *genjetWithNuSoftDrop_ = nullptr;
std::vector<const pat::Jet*> subjets_;
std::vector<const pat::PackedCandidate*> daughters_;
std::vector<const pat::Jet*> uncorr_subjets_;
std::vector<reco::CandidatePtr> daughters_;
std::map<reco::CandidatePtr::key_type, float> puppi_wgt_cache_;

};

Expand Down
54 changes: 35 additions & 19 deletions NtupleCommons/src/JetHelper.cc
View file
Expand Up @@ -9,46 +9,62 @@

namespace deepntuples {

JetHelper::JetHelper(const pat::Jet* jet, bool has_puppi_weighted_daughters) : jet_(jet), has_puppi_weighted_daughters_(has_puppi_weighted_daughters) {
JetHelper::JetHelper(const pat::Jet* jet, const edm::Handle<reco::CandidateView> &pfcands) : jet_(jet) {
if (!jet) throw cms::Exception("[JetHelper::JetHelper] Null pointer for input jet!");
if (jet->nSubjetCollections() == 0) throw cms::Exception("[JetHelper::JetHelper] No subjet collection for input jet!");
initializeConstituents();
initializeConstituents(pfcands);
}

void JetHelper::initializeConstituents() {
void JetHelper::initializeConstituents(const edm::Handle<reco::CandidateView> &pfcands) {
subjets_.clear();
uncorr_subjets_.clear();
daughters_.clear();
puppi_wgt_cache_.clear();

// get subjets
auto subjets = jet_->subjets();
for (const auto &sj : subjets){
subjets_.push_back(&(*sj));
uncorr_subjets_.push_back(&(*sj));
}
// sort subjets by pt
std::sort(subjets_.begin(), subjets_.end(),
[](const pat::Jet* p1, const pat::Jet* p2){return p1->pt()>p2->pt();});

std::sort(uncorr_subjets_.begin(), uncorr_subjets_.end(),
[](const pat::Jet* p1, const pat::Jet* p2){return p1->correctedP4("Uncorrected").pt()>p2->correctedP4("Uncorrected").pt();});

// get all consitituents
for (unsigned idau=0; idau<jet_->numberOfDaughters(); ++idau){
const auto *dau = jet_->daughter(idau);
if (dau->numberOfDaughters()>0){
// is a subjet; add all daughters
for (unsigned k=0; k<dau->numberOfDaughters(); ++k){
const auto *cand = dynamic_cast<const pat::PackedCandidate*>(dau->daughter(k));
auto dauPtr = jet_->daughterPtr(idau);
if (dauPtr->numberOfDaughters()>0){
// is a subjet
const auto *sj = dynamic_cast<const pat::Jet*>(&(*dauPtr));
if (!sj) throw cms::Exception("[JetHelper::initializeConstituents] Cannot convert to subjet!");
// add all daughters
for (unsigned k=0; k<sj->numberOfDaughters(); ++k){
const auto& candPtr = sj->daughterPtr(k);
const auto *cand = dynamic_cast<const pat::PackedCandidate*>(&(*candPtr));
if (cand->puppiWeight() < 0.01) continue; // [94X] ignore particles w/ extremely low puppi weights
daughters_.push_back(cand);
// Here we get the original PackedCandidate as stored in MiniAOD (i.e., not puppi weighted)
// https://github.com/cms-sw/cmssw/pull/28035
daughters_.push_back(pfcands->ptrAt(dauPtr.key()));
// For the Puppi weight, we get it from the new candidate in case it is recomputed
puppi_wgt_cache_[dauPtr.key()] = cand->puppiWeight();
}
}else{
const auto *cand = dynamic_cast<const pat::PackedCandidate*>(dau);
const auto& candPtr = dauPtr;
const auto *cand = dynamic_cast<const pat::PackedCandidate*>(&(*candPtr));
if (cand->puppiWeight() < 0.01) continue; // [94X] ignore particles w/ extremely low puppi weights
daughters_.push_back(cand);
// Here we get the original PackedCandidate as stored in MiniAOD (i.e., not puppi weighted)
// https://github.com/cms-sw/cmssw/pull/28035
daughters_.push_back(pfcands->ptrAt(dauPtr.key()));
// For the Puppi weight, we get it from the new candidate in case it is recomputed
puppi_wgt_cache_[dauPtr.key()] = cand->puppiWeight();
}
}
// sort by puppi weighted pt
if (has_puppi_weighted_daughters_){
std::sort(daughters_.begin(), daughters_.end(),
[](const pat::PackedCandidate* p1, const pat::PackedCandidate* p2){return p1->pt() > p2->pt();});
} else {
std::sort(daughters_.begin(), daughters_.end(),
[](const pat::PackedCandidate* p1, const pat::PackedCandidate* p2){return p1->puppiWeight()*p1->pt() > p2->puppiWeight()*p2->pt();});
}
// sort by original pt
std::sort(daughters_.begin(), daughters_.end(), [](const reco::CandidatePtr &a, const reco::CandidatePtr &b){ return a->pt() > b->pt(); });

}

Expand Down
9 changes: 6 additions & 3 deletions Ntupler/plugins/DeepNtuplizer.cc
View file
Expand Up @@ -38,9 +38,9 @@ class DeepNtuplizer : public edm::one::EDAnalyzer<edm::one::SharedResources> {
virtual void endJob() override;

double jetR = -1;
bool has_puppi_weighted_daughters_ = false;

edm::EDGetTokenT<edm::View<pat::Jet>> jetToken_;
edm::EDGetTokenT<edm::View<reco::Candidate>> candToken_;
edm::EDGetTokenT<edm::Association<reco::GenJetCollection>> genJetWithNuMatchToken_;
edm::EDGetTokenT<edm::Association<reco::GenJetCollection>> genJetWithNuSoftDropMatchToken_;

Expand All @@ -56,8 +56,8 @@ class DeepNtuplizer : public edm::one::EDAnalyzer<edm::one::SharedResources> {

DeepNtuplizer::DeepNtuplizer(const edm::ParameterSet& iConfig):
jetR(iConfig.getParameter<double>("jetR")),
has_puppi_weighted_daughters_(iConfig.getParameter<bool>("hasPuppiWeightedDaughters")),
jetToken_(consumes<edm::View<pat::Jet> >(iConfig.getParameter<edm::InputTag>("jets"))),
candToken_(consumes<edm::View<reco::Candidate>>(iConfig.getParameter<edm::InputTag>("pfcands"))),
genJetWithNuMatchToken_(consumes<edm::Association<reco::GenJetCollection>>(iConfig.getParameter<edm::InputTag>("genJetsMatch"))),
genJetWithNuSoftDropMatchToken_(consumes<edm::Association<reco::GenJetCollection>>(iConfig.getParameter<edm::InputTag>("genJetsSoftDropMatch")))
{
Expand Down Expand Up @@ -98,6 +98,9 @@ void DeepNtuplizer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSe
edm::Handle<edm::View<pat::Jet>> jets;
iEvent.getByToken(jetToken_, jets);

edm::Handle<edm::View<reco::Candidate>> candHandle;
iEvent.getByToken(candToken_, candHandle);

edm::Handle<edm::Association<reco::GenJetCollection>> genJetWithNuMatchHandle;
iEvent.getByToken(genJetWithNuMatchToken_, genJetWithNuMatchHandle);

Expand All @@ -108,7 +111,7 @@ void DeepNtuplizer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSe
bool write_ = true;

const auto& jet = jets->at(idx); // need to keep the JEC for puppi sdmass corr
JetHelper jet_helper(&jet, has_puppi_weighted_daughters_);
JetHelper jet_helper(&jet, candHandle);
jet_helper.setGenjetWithNu((*genJetWithNuMatchHandle)[jets->refAt(idx)]);
jet_helper.setGenjetWithNuSoftDrop((*genJetWithNuSoftDropMatchHandle)[jets->refAt(idx)]);

Expand Down
1 change: 1 addition & 0 deletions Ntupler/python/DeepNtuplizer_cfi.py
View file
@@ -1,6 +1,7 @@
import FWCore.ParameterSet.Config as cms

deepntuplizer = cms.EDAnalyzer('DeepNtuplizer',
pfcands = cms.InputTag("packedPFCandidates"),
vertices = cms.InputTag("offlineSlimmedPrimaryVertices"),
puInfo = cms.InputTag("slimmedAddPileupInfo"),
rhoInfo = cms.InputTag("fixedGridRhoFastjetAll"),
Expand Down
15 changes: 15 additions & 0 deletions Ntupler/run/2017/ak15/qcd-mg.conf
View file
@@ -0,0 +1,15 @@
/QCD_HT300to500_TuneCP5_13TeV-madgraph-pythia8/RunIIFall17MiniAODv2-PU2017_12Apr2018_94X_mc2017_realistic_v14-v1/MINIAODSIM
/QCD_HT500to700_TuneCP5_13TeV-madgraph-pythia8/RunIIFall17MiniAODv2-PU2017_12Apr2018_94X_mc2017_realistic_v14-v2/MINIAODSIM
/QCD_HT700to1000_TuneCP5_13TeV-madgraph-pythia8/RunIIFall17MiniAODv2-PU2017_12Apr2018_94X_mc2017_realistic_v14-v1/MINIAODSIM
/QCD_HT1000to1500_TuneCP5_13TeV-madgraph-pythia8/RunIIFall17MiniAODv2-PU2017_12Apr2018_94X_mc2017_realistic_v14-v1/MINIAODSIM
/QCD_HT1500to2000_TuneCP5_13TeV-madgraph-pythia8/RunIIFall17MiniAODv2-PU2017_12Apr2018_94X_mc2017_realistic_v14-v2/MINIAODSIM
/QCD_HT2000toInf_TuneCP5_13TeV-madgraph-pythia8/RunIIFall17MiniAODv2-PU2017_12Apr2018_94X_mc2017_realistic_v14-v2/MINIAODSIM

#/QCD_HT100to200_TuneCP5_13TeV-madgraph-pythia8/RunIIFall17MiniAODv2-PU2017_12Apr2018_94X_mc2017_realistic_v14_ext1-v1/MINIAODSIM
#/QCD_HT200to300_TuneCP5_13TeV-madgraph-pythia8/RunIIFall17MiniAODv2-PU2017_12Apr2018_new_pmx_94X_mc2017_realistic_v14-v2/MINIAODSIM
#/QCD_HT300to500_TuneCP5_13TeV-madgraph-pythia8/RunIIFall17MiniAODv2-PU2017_12Apr2018_new_pmx_94X_mc2017_realistic_v14-v1/MINIAODSIM
#/QCD_HT500to700_TuneCP5_13TeV-madgraph-pythia8/RunIIFall17MiniAODv2-PU2017_12Apr2018_old_pmx_94X_mc2017_realistic_v14-v1/MINIAODSIM
#/QCD_HT700to1000_TuneCP5_13TeV-madgraph-pythia8/RunIIFall17MiniAODv2-PU2017_12Apr2018_new_pmx_94X_mc2017_realistic_v14-v2/MINIAODSIM
#/QCD_HT1000to1500_TuneCP5_13TeV-madgraph-pythia8/RunIIFall17MiniAODv2-PU2017_12Apr2018_new_pmx_94X_mc2017_realistic_v14-v1/MINIAODSIM
#/QCD_HT1500to2000_TuneCP5_13TeV-madgraph-pythia8/RunIIFall17MiniAODv2-PU2017_12Apr2018_old_pmx_94X_mc2017_realistic_v14-v1/MINIAODSIM
#/QCD_HT2000toInf_TuneCP5_13TeV-madgraph-pythia8/RunIIFall17MiniAODv2-PU2017_12Apr2018_old_pmx_94X_mc2017_realistic_v14-v1/MINIAODSIM

0 comments on commit bb65a79

Please sign in to comment.