CMS 3D CMS Logo

VBFGenJetFilter.cc
Go to the documentation of this file.
2 
5 
6 #include <HepMC/GenVertex.h>
7 
8 // ROOT includes
9 #include "TMath.h"
10 
11 // C++ includes
12 #include <iostream>
13 
14 using namespace edm;
15 using namespace std;
16 
17 
19 oppositeHemisphere(iConfig.getUntrackedParameter<bool> ("oppositeHemisphere",false)),
20 leadJetsNoLepMass (iConfig.getUntrackedParameter<bool> ("leadJetsNoLepMass", false)),
21 ptMin (iConfig.getUntrackedParameter<double>("minPt", 20)),
22 etaMin (iConfig.getUntrackedParameter<double>("minEta", -5.0)),
23 etaMax (iConfig.getUntrackedParameter<double>("maxEta", 5.0)),
24 minInvMass (iConfig.getUntrackedParameter<double>("minInvMass", 0.0)),
25 maxInvMass (iConfig.getUntrackedParameter<double>("maxInvMass", 99999.0)),
26 minDeltaPhi (iConfig.getUntrackedParameter<double>("minDeltaPhi", -1.0)),
27 maxDeltaPhi (iConfig.getUntrackedParameter<double>("maxDeltaPhi", 99999.0)),
28 minDeltaEta (iConfig.getUntrackedParameter<double>("minDeltaEta", -1.0)),
29 maxDeltaEta (iConfig.getUntrackedParameter<double>("maxDeltaEta", 99999.0)),
30 minLeadingJetsInvMass (iConfig.getUntrackedParameter<double>("minLeadingJetsInvMass", 0.0)),
31 maxLeadingJetsInvMass (iConfig.getUntrackedParameter<double>("maxLeadingJetsInvMass", 99999.0)),
32 deltaRJetLep (iConfig.getUntrackedParameter<double>("deltaRJetLep", 0.3))
33 {
34 
35  m_inputTag_GenJetCollection = consumes<reco::GenJetCollection>(iConfig.getUntrackedParameter<edm::InputTag>("inputTag_GenJetCollection",edm::InputTag("ak5GenJetsNoNu")));
36  if (leadJetsNoLepMass) m_inputTag_GenParticleCollection = consumes<reco::GenParticleCollection>(iConfig.getUntrackedParameter<edm::InputTag>("genParticles",edm::InputTag("genParticles")));
37 
38 }
39 
41 
42 }
43 
44 
45 vector<const reco::GenParticle*> VBFGenJetFilter::filterGenLeptons(const vector<reco::GenParticle>* particles){
46  vector<const reco::GenParticle*> out;
47 
48 
49 
50  for(const auto & p : *particles){
51 
52  int absPdgId = std::abs(p.pdgId());
53 
54  if(((absPdgId == 11) || (absPdgId == 13) || (absPdgId == 15)) && p.isHardProcess()) {
55  out.push_back(&p);
56  }
57 
58 
59  }
60  return out;
61 }
62 
63 
64 
65 vector<const reco::GenJet*> VBFGenJetFilter::filterGenJets(const vector<reco::GenJet>* jets){
66 
67  vector<const reco::GenJet*> out;
68 
69  for(unsigned i=0; i<jets->size(); i++){
70 
71  const reco::GenJet* j = &((*jets)[i]);
72 
73  if(j->p4().pt() >ptMin && j->p4().eta()>etaMin && j->p4().eta()<etaMax)
74  {
75  out.push_back(j);
76  }
77  }
78 
79  return out;
80 }
81 
82 
83 // ------------ method called to skim the data ------------
85 {
86  using namespace edm;
87 
88  Handle< vector<reco::GenJet> > handleGenJets;
89  iEvent.getByToken(m_inputTag_GenJetCollection, handleGenJets);
90  const vector<reco::GenJet>* genJets = handleGenJets.product();
91 
92  // Getting filtered generator jets
93  vector<const reco::GenJet*> filGenJets = filterGenJets(genJets);
94 
95  // If we do not find at least 2 jets veto the event
96  if(filGenJets.size()<2){return false;}
97 
98 
99  // Testing dijet mass
100  if(leadJetsNoLepMass) {
101 
102  Handle<reco::GenParticleCollection> genParticelesCollection;
103  iEvent.getByToken(m_inputTag_GenParticleCollection, genParticelesCollection);
104  const vector<reco::GenParticle>* genParticles = genParticelesCollection.product();
105 
106  // Getting filtered generator muons
107  vector<const reco::GenParticle*> filGenLep = filterGenLeptons(genParticles);
108 
109  // Getting p4 of jet with no lepton
110  vector<math::XYZTLorentzVector> genJetsWithoutLeptonsP4;
111  unsigned int jetIdx = 0;
112 
113  while(genJetsWithoutLeptonsP4.size()<2 && jetIdx < filGenJets.size()) {
114  bool jetWhitoutLep = true;
115 
116  const math::XYZTLorentzVector & p4J= (filGenJets[jetIdx])->p4();
117  for(unsigned int i = 0; i < filGenLep.size() && jetWhitoutLep; ++i) {
118  if(reco::deltaR2((filGenLep[i])->p4(), p4J) < deltaRJetLep*deltaRJetLep)
119  jetWhitoutLep = false;
120  }
121  if (jetWhitoutLep) genJetsWithoutLeptonsP4.push_back(p4J);
122  ++jetIdx;
123  }
124 
125  // Checking the invariant mass of the leading jets
126  if (genJetsWithoutLeptonsP4.size() < 2) return false;
127  float invMassLeadingJet = (genJetsWithoutLeptonsP4[0] + genJetsWithoutLeptonsP4[1]).M();
128  if ( invMassLeadingJet > minLeadingJetsInvMass && invMassLeadingJet < maxLeadingJetsInvMass) return true;
129  else return false;
130  }
131 
132 
133  for(unsigned a=0; a<filGenJets.size(); a++){
134  for(unsigned b=a+1; b<filGenJets.size(); b++){
135 
136  const reco::GenJet* pA = filGenJets[a];
137  const reco::GenJet* pB = filGenJets[b];
138 
139  // Getting the dijet vector
140  math::XYZTLorentzVector diJet = pA->p4() + pB->p4();
141 
142  // Testing opposite hemispheres
143  double dijetProd = pA->p4().eta()*pB->p4().eta();
144  if(oppositeHemisphere && dijetProd>=0){continue;}
145 
146  // Testing dijet mass
147  double invMass = diJet.mass();
148  if(invMass<=minInvMass || invMass>maxInvMass){continue;}
149 
150  // Testing dijet delta eta
151  double dEta = fabs(pA->p4().eta()-pB->p4().eta());
152  if(dEta<=minDeltaEta || dEta>maxDeltaEta){continue;}
153 
154  // Testing dijet delta phi
155  double dPhi = fabs(reco::deltaPhi(pA->p4().phi(),pB->p4().phi()));
156  if(dPhi<=minDeltaPhi || dPhi>maxDeltaPhi){continue;}
157 
158  return true;
159  }
160  }
161 
162  return false;
163 }
164 
165 //define this as a plug-in
constexpr double deltaPhi(double phi1, double phi2)
Definition: deltaPhi.h:22
T getUntrackedParameter(std::string const &, T const &) const
double minLeadingJetsInvMass
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
std::vector< const reco::GenJet * > filterGenJets(const std::vector< reco::GenJet > *jets)
bool filter(edm::Event &, const edm::EventSetup &) override
edm::EDGetTokenT< reco::GenParticleCollection > m_inputTag_GenParticleCollection
VBFGenJetFilter(const edm::ParameterSet &)
std::vector< const reco::GenParticle * > filterGenLeptons(const std::vector< reco::GenParticle > *particles)
double maxLeadingJetsInvMass
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
int iEvent
Definition: GenABIO.cc:224
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
double p4[4]
Definition: TauolaWrapper.h:92
vector< PseudoJet > jets
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
Jets made from MC generator particles.
Definition: GenJet.h:25
const LorentzVector & p4() const final
four-momentum Lorentz vector
Definition: LeafCandidate.h:99
constexpr auto deltaR2(const T1 &t1, const T2 &t2) -> decltype(t1.eta())
Definition: deltaR.h:16
T const * product() const
Definition: Handle.h:74
double b
Definition: hdecay.h:120
HLT enums.
double a
Definition: hdecay.h:121
~VBFGenJetFilter() override
edm::EDGetTokenT< reco::GenJetCollection > m_inputTag_GenJetCollection