CMS 3D CMS Logo

HLTVertexFilter.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: HLTVertexFilter
4 // Class: HLTVertexFilter
5 //
14 // Original Author: Andrea Bocci
15 // Created: Tue Apr 20 12:34:27 CEST 2010
16 
17 
18 #include <cmath>
19 
20 // user include files
32 
33 //
34 // class declaration
35 //
36 class HLTVertexFilter : public HLTFilter {
37 public:
38  explicit HLTVertexFilter(const edm::ParameterSet & config);
39  ~HLTVertexFilter() override;
40  static void fillDescriptions(edm::ConfigurationDescriptions & descriptions);
41 
42 private:
43 
44  bool hltFilter(edm::Event & event, const edm::EventSetup & setup, trigger::TriggerFilterObjectWithRefs & filterproduct) const override;
45 
47  edm::InputTag m_inputTag; // input vertex collection
48  double m_minNDoF; // minimum vertex NDoF
49  double m_maxChi2; // maximum vertex chi2
50  double m_maxD0; // maximum transverse distance from the beam
51  double m_maxZ; // maximum longitudinal distance nominal center of the detector
52  unsigned int m_minVertices;
53 
54 };
55 
56 //
57 // constructors and destructor
58 //
60  m_inputTag(config.getParameter<edm::InputTag>("inputTag")),
61  m_minNDoF(config.getParameter<double>("minNDoF")),
62  m_maxChi2(config.getParameter<double>("maxChi2")),
63  m_maxD0(config.getParameter<double>("maxD0")),
64  m_maxZ(config.getParameter<double>("maxZ")),
65  m_minVertices(config.getParameter<unsigned int>("minVertices"))
66 {
67  m_inputToken = consumes<reco::VertexCollection>(m_inputTag);
68 }
69 
70 
72 
73 void
77  desc.add<edm::InputTag>("inputTag",edm::InputTag("hltPixelVertices"));
78  desc.add<double>("minNDoF",0.);
79  desc.add<double>("maxChi2",99999.);
80  desc.add<double>("maxD0",1.);
81  desc.add<double>("maxZ",15.);
82  desc.add<unsigned int>("minVertices",1);
83  descriptions.add("hltVertexFilter",desc);
84 }
85 
86 
87 //
88 // member functions
89 //
90 
91 // ------------ method called on each new Event ------------
92 bool
94 
95  // get hold of collection of objects
97  event.getByToken(m_inputToken, vertices);
98 
99  unsigned int n = 0;
100  if (vertices.isValid()) {
101  for(auto const& vertex : * vertices) {
102  if (vertex.isValid()
103  and not vertex.isFake()
104  and (vertex.chi2() <= m_maxChi2)
105  and (vertex.ndof() >= m_minNDoF)
106  and (std::abs(vertex.z()) <= m_maxZ)
107  and (std::abs(vertex.y()) <= m_maxD0)
108  and (std::abs(vertex.x()) <= m_maxD0)
109  and (std::sqrt(std::pow(vertex.x(),2)+std::pow(vertex.y(),2)) <= m_maxD0)
110  )
111  ++n;
112  }
113  }
114 
115  // filter decision
116  return (n >= m_minVertices);
117 }
118 
119 // declare this class as a framework plugin
edm::InputTag m_inputTag
def setup(process, global_tag, zero_tesla=False)
Definition: GeneralSetup.py:2
Definition: config.py:1
~HLTVertexFilter() override
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
T sqrt(T t)
Definition: SSEVec.h:18
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
ParameterDescriptionBase * add(U const &iLabel, T const &value)
edm::EDGetTokenT< reco::VertexCollection > m_inputToken
bool isValid() const
Definition: HandleBase.h:74
static void makeHLTFilterDescription(edm::ParameterSetDescription &desc)
Definition: HLTFilter.cc:29
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
void add(std::string const &label, ParameterSetDescription const &psetDescription)
unsigned int m_minVertices
HLT enums.
HLTVertexFilter(const edm::ParameterSet &config)
bool hltFilter(edm::Event &event, const edm::EventSetup &setup, trigger::TriggerFilterObjectWithRefs &filterproduct) const override
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
Definition: event.py:1