diff --git a/GeneratorInterface/Core/interface/EmbeddingHepMCFilter.h b/GeneratorInterface/Core/interface/EmbeddingHepMCFilter.h index 2ce970a6e1ea1..17f7a4e6fa29a 100644 --- a/GeneratorInterface/Core/interface/EmbeddingHepMCFilter.h +++ b/GeneratorInterface/Core/interface/EmbeddingHepMCFilter.h @@ -4,6 +4,7 @@ #include "FWCore/MessageLogger/interface/MessageLogger.h" #include "GeneratorInterface/Core/interface/BaseHepMCFilter.h" +#include "PhysicsTools/HepMCCandAlgos/interface/MCTruthHelper.h" #include "DataFormats/Candidate/interface/Candidate.h" class EmbeddingHepMCFilter : public BaseHepMCFilter { @@ -15,6 +16,8 @@ class EmbeddingHepMCFilter : public BaseHepMCFilter { const int electron_neutrino_PDGID_ = 12; const int electronPDGID_ = 11; int ZPDGID_ = 23; + bool includeDY_ = false; + MCTruthHelper mcTruthHelper_; enum class TauDecayMode : int { Unfilled = -1, Muon = 0, Electron = 1, Hadronic = 2 }; diff --git a/GeneratorInterface/Core/src/EmbeddingHepMCFilter.cc b/GeneratorInterface/Core/src/EmbeddingHepMCFilter.cc index 499ae555ac108..1d9a21b9eaf08 100644 --- a/GeneratorInterface/Core/src/EmbeddingHepMCFilter.cc +++ b/GeneratorInterface/Core/src/EmbeddingHepMCFilter.cc @@ -4,7 +4,8 @@ #include "boost/algorithm/string/trim_all.hpp" EmbeddingHepMCFilter::EmbeddingHepMCFilter(const edm::ParameterSet &iConfig) - : ZPDGID_(iConfig.getParameter("BosonPDGID")) { + : ZPDGID_(iConfig.getParameter("BosonPDGID")), + includeDY_(iConfig.existsAs("IncludeDY") ? iConfig.getParameter("IncludeDY") : false) { // Defining standard decay channels ee.fill(TauDecayMode::Electron); ee.fill(TauDecayMode::Electron); @@ -63,25 +64,35 @@ bool EmbeddingHepMCFilter::filter(const HepMC::GenEvent *evt) { // One can stop the loop after the second tau is reached and processed. for (HepMC::GenEvent::particle_const_iterator particle = evt->particles_begin(); particle != evt->particles_end(); ++particle) { - int mom_id = 0; // No particle available with PDG ID 0 - if ((*particle)->production_vertex() != nullptr) { // search for the mom via the production_vertex - if ((*particle)->production_vertex()->particles_in_const_begin() != - (*particle)->production_vertex()->particles_in_const_end()) { - mom_id = (*(*particle)->production_vertex()->particles_in_const_begin())->pdg_id(); // mom was found + int mom_id = 0; // no particle available with PDG ID 0 + bool isHardProc = false; // mother is ZPDGID_, or is lepton from hard process (DY process qq -> ll) + int pdg_id = std::abs((*particle)->pdg_id()); + HepMC::GenVertex *vertex = (*particle)->production_vertex(); + if (vertex != nullptr) { // search for the mom via the production_vertex + HepMC::GenVertex::particles_in_const_iterator mom = vertex->particles_in_const_begin(); + if (mom != vertex->particles_in_const_end()) { + mom_id = std::abs((*mom)->pdg_id()); // mom was found + } + if (mom_id == ZPDGID_) { + isHardProc = true; // intermediate boson + } else if (includeDY_ && 11 <= pdg_id && pdg_id <= 16 && mcTruthHelper_.isFirstCopy(**particle) && + mcTruthHelper_.fromHardProcess(**particle)) { + edm::LogInfo("EmbeddingHepMCFilter") << (*particle)->pdg_id() << " with mother " << (*mom)->pdg_id(); + isHardProc = true; // assume Drell-Yan qq -> ll without intermediate boson } } - if (std::abs((*particle)->pdg_id()) == tauonPDGID_ && mom_id == ZPDGID_) { + if (!isHardProc) { + continue; + } else if (pdg_id == tauonPDGID_) { reco::Candidate::LorentzVector p4Vis; decay_and_sump4Vis((*particle), p4Vis); // recursive access to final states. p4VisPair_.push_back(p4Vis); - } else if (std::abs((*particle)->pdg_id()) == muonPDGID_ && - mom_id == ZPDGID_) { // Also handle the option when Z-> mumu + } else if (pdg_id == muonPDGID_) { // Also handle the option when Z-> mumu reco::Candidate::LorentzVector p4Vis = (reco::Candidate::LorentzVector)(*particle)->momentum(); DecayChannel_.fill(TauDecayMode::Muon); // take the muon cuts p4VisPair_.push_back(p4Vis); - } else if (std::abs((*particle)->pdg_id()) == electronPDGID_ && - mom_id == ZPDGID_) { // Also handle the option when Z-> ee + } else if (pdg_id == electronPDGID_) { // Also handle the option when Z-> ee reco::Candidate::LorentzVector p4Vis = (reco::Candidate::LorentzVector)(*particle)->momentum(); DecayChannel_.fill(TauDecayMode::Electron); // take the electron cuts p4VisPair_.push_back(p4Vis); diff --git a/TauAnalysis/MCEmbeddingTools/python/EmbeddingPythia8Hadronizer_cfi.py b/TauAnalysis/MCEmbeddingTools/python/EmbeddingPythia8Hadronizer_cfi.py index 5c056288c8203..ef64aad95b4de 100644 --- a/TauAnalysis/MCEmbeddingTools/python/EmbeddingPythia8Hadronizer_cfi.py +++ b/TauAnalysis/MCEmbeddingTools/python/EmbeddingPythia8Hadronizer_cfi.py @@ -19,7 +19,8 @@ MuHadCut = cms.string('Mu.Pt > 18 && Had.Pt > 25 && Mu.Eta < 2.1'), MuMuCut = cms.string('Mu1.Pt > 17 && Mu2.Pt > 8'), Final_States = cms.vstring('ElEl','ElHad','ElMu','HadHad','MuHad','MuMu'), - BosonPDGID = cms.int32(23) + BosonPDGID = cms.int32(23), + IncludeDY = cms.bool(False) ), ), pythiaPylistVerbosity = cms.untracked.int32(0),