CMS 3D CMS Logo

HLTJetVBFFilter.cc
Go to the documentation of this file.
1 
20 
21 //
22 // constructors and destructor
23 //
24 template<typename T>
26  inputTag_ (iConfig.template getParameter< edm::InputTag > ("inputTag")),
27  minPtLow_ (iConfig.template getParameter<double> ("minPtLow")),
28  minPtHigh_ (iConfig.template getParameter<double> ("minPtHigh")),
29  etaOpposite_ (iConfig.template getParameter<bool> ("etaOpposite")),
30  minDeltaEta_ (iConfig.template getParameter<double> ("minDeltaEta")),
31  minInvMass_ (iConfig.template getParameter<double> ("minInvMass")),
32  maxEta_ (iConfig.template getParameter<double> ("maxEta")),
33  leadingJetOnly_ (iConfig.template getParameter<bool> ("leadingJetOnly")),
34  triggerType_ (iConfig.template getParameter<int> ("triggerType"))
35 {
36  m_theObjectToken = consumes<std::vector<T>>(inputTag_);
37  LogDebug("") << "HLTJetVBFFilter: Input/minPtLow_/minPtHigh_/triggerType : "
38  << inputTag_.encode() << " "
39  << minPtLow_ << " "
40  << minPtHigh_ << " "
41  << triggerType_;
42 }
43 
44 template<typename T>
46 
47 template<typename T>
48 void
52  desc.add<edm::InputTag>("inputTag",edm::InputTag("hltAntiKT5ConvPFJets"));
53  desc.add<double>("minPtLow",40.);
54  desc.add<double>("minPtHigh",40.);
55  desc.add<bool>("etaOpposite",false);
56  desc.add<double>("minDeltaEta",4.0);
57  desc.add<double>("minInvMass",1000.);
58  desc.add<double>("maxEta",5.0);
59  desc.add<bool>("leadingJetOnly",false);
60  desc.add<int>("triggerType",trigger::TriggerJet);
61  descriptions.add(defaultModuleLabel<HLTJetVBFFilter<T>>(), desc);
62 }
63 
64 //
65 // ------------ method called to produce the data ------------
66 //
67 template<typename T>
68 bool
70 {
71  using namespace std;
72  using namespace edm;
73  using namespace reco;
74  using namespace trigger;
75 
76  typedef vector<T> TCollection;
77  typedef Ref<TCollection> TRef;
78 
79  // The filter object
80  if (saveTags()) filterproduct.addCollectionTag(inputTag_);
81 
82  // get hold of collection of objects
84  iEvent.getByToken (m_theObjectToken,objects);
85 
86  // look at all candidates, check cuts and add to filter object
87  int n(0);
88 
89  // events with two or more jets
90  if(objects->size() > 1){
91 
92  double ejet1 = 0.;
93  double pxjet1 = 0.;
94  double pyjet1 = 0.;
95  double pzjet1 = 0.;
96  double ptjet1 = 0.;
97  double etajet1 = 0.;
98 
99  double ejet2 = 0.;
100  double pxjet2 = 0.;
101  double pyjet2 = 0.;
102  double pzjet2 = 0.;
103  double ptjet2 = 0.;
104  double etajet2 = 0.;
105 
106  // loop on all jets
107  int countJet1(0);
108  int countJet2(0);
109  typename TCollection::const_iterator jet1 ( objects->begin() );
110  for (; jet1!=objects->end(); jet1++) {
111  countJet1++;
112  if( leadingJetOnly_==true && countJet1>2 ) break;
113  //
114  if( jet1->pt() < minPtHigh_ ) break; //No need to go to the next jet (lower PT)
115  if( std::abs(jet1->eta()) > maxEta_ ) continue;
116  //
117  countJet2 = countJet1-1;
118  typename TCollection::const_iterator jet2 ( jet1+1 );
119  for (; jet2!=objects->end(); jet2++) {
120  countJet2++;
121  if( leadingJetOnly_==true && countJet2>2 ) break;
122  //
123  if( jet2->pt() < minPtLow_ ) break; //No need to go to the next jet (lower PT)
124  if( std::abs(jet2->eta()) > maxEta_ ) continue;
125  //
126  ejet1 = jet1->energy();
127  pxjet1 = jet1->px();
128  pyjet1 = jet1->py();
129  pzjet1 = jet1->pz();
130  ptjet1 = jet1->pt();
131  etajet1 = jet1->eta();
132 
133  ejet2 = jet2->energy();
134  pxjet2 = jet2->px();
135  pyjet2 = jet2->py();
136  pzjet2 = jet2->pz();
137  ptjet2 = jet2->pt();
138  etajet2 = jet2->eta();
139  //
140  float deltaetajet = etajet1 - etajet2;
141  float invmassjet = sqrt( (ejet1 + ejet2) * (ejet1 + ejet2) -
142  (pxjet1 + pxjet2) * (pxjet1 + pxjet2) -
143  (pyjet1 + pyjet2) * (pyjet1 + pyjet2) -
144  (pzjet1 + pzjet2) * (pzjet1 + pzjet2) );
145 
146  // VBF cuts
147  if ( (ptjet1 > minPtHigh_) &&
148  (ptjet2 > minPtLow_) &&
149  ( (etaOpposite_ == true && etajet1*etajet2 < 0) || (etaOpposite_ == false) ) &&
150  (std::abs(deltaetajet) > minDeltaEta_) &&
151  (std::abs(invmassjet) > minInvMass_) ){
152  ++n;
153  TRef ref1 = TRef(objects,distance(objects->begin(),jet1));
154  TRef ref2 = TRef(objects,distance(objects->begin(),jet2));
155  filterproduct.addObject(triggerType_,ref1);
156  filterproduct.addObject(triggerType_,ref2);
157  }// VBF cuts
158  //if(n>=1) break; //Store all possible pairs
159  }
160  //if(n>=1) break; //Store all possible pairs
161  }// loop on all jets
162  }// events with two or more jets
163 
164  // filter decision
165  bool accept(n>=1);
166 
167  return accept;
168 }
#define LogDebug(id)
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
std::string defaultModuleLabel()
bool accept(const edm::Event &event, const edm::TriggerResults &triggerTable, const std::string &triggerPath)
Definition: TopDQMHelpers.h:30
std::string encode() const
Definition: InputTag.cc:159
void addObject(int id, const reco::RecoEcalCandidateRef &ref)
setters for L3 collections: (id=physics type, and Ref<C>)
int iEvent
Definition: GenABIO.cc:224
edm::InputTag inputTag_
T sqrt(T t)
Definition: SSEVec.h:18
~HLTJetVBFFilter() override
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
def template(fileName, svg, replaceme="REPLACEME")
Definition: svgfig.py:521
ParameterDescriptionBase * add(U const &iLabel, T const &value)
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
static void makeHLTFilterDescription(edm::ParameterSetDescription &desc)
Definition: HLTFilter.cc:29
HLTJetVBFFilter(const edm::ParameterSet &)
void addCollectionTag(const edm::InputTag &collectionTag)
collectionTags
void add(std::string const &label, ParameterSetDescription const &psetDescription)
fixed size matrix
bool saveTags() const
Definition: HLTFilter.h:45
HLT enums.
edm::EDGetTokenT< std::vector< T > > m_theObjectToken
bool hltFilter(edm::Event &, const edm::EventSetup &, trigger::TriggerFilterObjectWithRefs &filterproduct) const override