File indexing completed on 2024-04-06 12:18:43
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017 #include <cmath>
0018
0019
0020 #include "FWCore/Framework/interface/Frameworkfwd.h"
0021 #include "FWCore/Framework/interface/Event.h"
0022 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0023 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0024 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0025 #include "DataFormats/Common/interface/Handle.h"
0026 #include "DataFormats/VertexReco/interface/Vertex.h"
0027 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0028 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0029 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0030 #include "HLTrigger/HLTcore/interface/HLTFilter.h"
0031
0032
0033
0034
0035 class HLTVertexFilter : public HLTFilter {
0036 public:
0037 explicit HLTVertexFilter(const edm::ParameterSet& config);
0038 ~HLTVertexFilter() override;
0039 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0040
0041 private:
0042 bool hltFilter(edm::Event& event,
0043 const edm::EventSetup& setup,
0044 trigger::TriggerFilterObjectWithRefs& filterproduct) const override;
0045
0046 edm::EDGetTokenT<reco::VertexCollection> m_inputToken;
0047 edm::InputTag m_inputTag;
0048 double m_minNDoF;
0049 double m_maxChi2;
0050 double m_maxD0;
0051 double m_maxZ;
0052 unsigned int m_minVertices;
0053 };
0054
0055
0056
0057
0058 HLTVertexFilter::HLTVertexFilter(const edm::ParameterSet& config)
0059 : HLTFilter(config),
0060 m_inputTag(config.getParameter<edm::InputTag>("inputTag")),
0061 m_minNDoF(config.getParameter<double>("minNDoF")),
0062 m_maxChi2(config.getParameter<double>("maxChi2")),
0063 m_maxD0(config.getParameter<double>("maxD0")),
0064 m_maxZ(config.getParameter<double>("maxZ")),
0065 m_minVertices(config.getParameter<unsigned int>("minVertices")) {
0066 m_inputToken = consumes<reco::VertexCollection>(m_inputTag);
0067 }
0068
0069 HLTVertexFilter::~HLTVertexFilter() = default;
0070
0071 void HLTVertexFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0072 edm::ParameterSetDescription desc;
0073 makeHLTFilterDescription(desc);
0074 desc.add<edm::InputTag>("inputTag", edm::InputTag("hltPixelVertices"));
0075 desc.add<double>("minNDoF", 0.);
0076 desc.add<double>("maxChi2", 99999.);
0077 desc.add<double>("maxD0", 1.);
0078 desc.add<double>("maxZ", 15.);
0079 desc.add<unsigned int>("minVertices", 1);
0080 descriptions.add("hltVertexFilter", desc);
0081 }
0082
0083
0084
0085
0086
0087
0088 bool HLTVertexFilter::hltFilter(edm::Event& event,
0089 edm::EventSetup const& setup,
0090 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0091
0092 edm::Handle<reco::VertexCollection> vertices;
0093 event.getByToken(m_inputToken, vertices);
0094
0095 unsigned int n = 0;
0096 if (vertices.isValid()) {
0097 for (auto const& vertex : *vertices) {
0098 if (vertex.isValid() and not vertex.isFake() and (vertex.chi2() <= m_maxChi2) and (vertex.ndof() >= m_minNDoF) and
0099 (std::abs(vertex.z()) <= m_maxZ) and (std::abs(vertex.y()) <= m_maxD0) and
0100 (std::abs(vertex.x()) <= m_maxD0) and
0101 (std::sqrt(std::pow(vertex.x(), 2) + std::pow(vertex.y(), 2)) <= m_maxD0))
0102 ++n;
0103 }
0104 }
0105
0106
0107 return (n >= m_minVertices);
0108 }
0109
0110
0111 #include "FWCore/Framework/interface/MakerMacros.h"
0112 DEFINE_FWK_MODULE(HLTVertexFilter);