CMS 3D CMS Logo

HLTJetPairDzMatchFilter.cc
Go to the documentation of this file.
1 
15 
16 template<typename T>
18 {
19  m_jetTag = conf.getParameter<edm::InputTag>("JetSrc");
20  m_jetToken = consumes<std::vector<T> >(m_jetTag);
21  m_jetMinPt = conf.getParameter<double>("JetMinPt");
22  m_jetMaxEta = conf.getParameter<double>("JetMaxEta");
23  m_jetMinDR = conf.getParameter<double>("JetMinDR");
24  m_jetMaxDZ = conf.getParameter<double>("JetMaxDZ");
25  m_triggerType = conf.getParameter<int>("TriggerType");
26 
27  // set the minimum DR between jets, so that one never has a chance
28  // to create a pair out of the same Jet replicated with two
29  // different vertices
30  if (m_jetMinDR < 0.1) m_jetMinDR = 0.1;
31 
32 }
33 
34 template<typename T>
35 void
39  desc.add<edm::InputTag>("JetSrc",edm::InputTag("hltMatchL2Tau30ToPixelTrk5"));
40  desc.add<double>("JetMinPt",25.0);
41  desc.add<double>("JetMaxEta",2.4);
42  desc.add<double>("JetMinDR",0.2);
43  desc.add<double>("JetMaxDZ",0.2);
44  desc.add<int>("TriggerType",84);
45  descriptions.add(defaultModuleLabel<HLTJetPairDzMatchFilter<T>>(), desc);
46 }
47 
48 template<typename T>
50 
51 template<typename T>
53 {
54  using namespace std;
55  using namespace edm;
56  using namespace reco;
57 
58  typedef vector<T> TCollection;
59  typedef Ref<TCollection> TRef;
60 
61  // The resuilting filter object to store in the Event
62 
63  if (saveTags()) filterproduct.addCollectionTag(m_jetTag);
64 
65  // Ref to Candidate object to be recorded in the filter object
66  TRef ref;
67 
68  // *** Pick up L2 tau jets which have been equipped with some meaningful vertices before that ***
69 
70  edm::Handle<TCollection> jetsHandle;
71  ev.getByToken(m_jetToken, jetsHandle );
72  const TCollection & jets = *jetsHandle;
73  const size_t n_jets = jets.size();
74 
75  // *** Combine jets into pairs and check the dz matching ***
76 
77  size_t npairs = 0, nfail_dz = 0;
78  if (n_jets > 1) for(size_t j1 = 0; j1 < n_jets; ++j1)
79  {
80  if ( jets[j1].pt() < m_jetMinPt || std::abs(jets[j1].eta()) > m_jetMaxEta ) continue;
81 
82  float mindz = 99.f;
83  for(size_t j2 = j1+1; j2 < n_jets; ++j2)
84  {
85  if ( jets[j2].pt() < m_jetMinPt || std::abs(jets[j2].eta()) > m_jetMaxEta ) continue;
86 
87  float deta = jets[j1].eta() - jets[j2].eta();
88  float dphi = reco::deltaPhi(jets[j1].phi(), jets[j2].phi());
89  float dr2 = dphi*dphi + deta*deta;
90  float dz = jets[j1].vz() - jets[j2].vz();
91 
92  // skip pairs of jets that are close
93  if ( dr2 < m_jetMinDR*m_jetMinDR ) {
94  continue;
95  }
96 
97  if (std::abs(dz) < std::abs(mindz)) mindz = dz;
98 
99  // do not form a pair if dz is too large
100  if ( std::abs(dz) > m_jetMaxDZ ) {
101  ++nfail_dz;
102  continue;
103  }
104 
105  // add references to both jets
106  ref = TRef(jetsHandle, j1);
107  filterproduct.addObject(m_triggerType, ref);
108 
109  ref = TRef(jetsHandle, j2);
110  filterproduct.addObject(m_triggerType, ref);
111 
112  ++npairs;
113 
114  }
115 
116  }
117 
118  return npairs>0;
119 
120 }
constexpr double deltaPhi(double phi1, double phi2)
Definition: deltaPhi.h:22
T getParameter(std::string const &) const
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
bool ev
std::string defaultModuleLabel()
void addObject(int id, const reco::RecoEcalCandidateRef &ref)
setters for L3 collections: (id=physics type, and Ref<C>)
vector< PseudoJet > jets
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
bool hltFilter(edm::Event &, const edm::EventSetup &, trigger::TriggerFilterObjectWithRefs &filterproduct) const override
ParameterDescriptionBase * add(U const &iLabel, T const &value)
static void makeHLTFilterDescription(edm::ParameterSetDescription &desc)
Definition: HLTFilter.cc:29
void addCollectionTag(const edm::InputTag &collectionTag)
collectionTags
void add(std::string const &label, ParameterSetDescription const &psetDescription)
HLTJetPairDzMatchFilter(const edm::ParameterSet &)
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
fixed size matrix
bool saveTags() const
Definition: HLTFilter.h:45
HLT enums.
edm::EDGetTokenT< std::vector< T > > m_jetToken
~HLTJetPairDzMatchFilter() override