File indexing completed on 2023-03-17 10:53:45
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021 #include <memory>
0022
0023
0024 #include "FWCore/Framework/interface/Frameworkfwd.h"
0025 #include "FWCore/Framework/interface/global/EDFilter.h"
0026
0027 #include "FWCore/Framework/interface/Event.h"
0028 #include "FWCore/Framework/interface/MakerMacros.h"
0029
0030 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0031 #include "FWCore/Utilities/interface/InputTag.h"
0032 #include "FWCore/Utilities/interface/EDGetToken.h"
0033 #include "DataFormats/VertexReco/interface/Vertex.h"
0034 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0035
0036
0037
0038
0039 class GoodVertexFilter : public edm::global::EDFilter<> {
0040 public:
0041 explicit GoodVertexFilter(const edm::ParameterSet&);
0042 ~GoodVertexFilter() override;
0043
0044 private:
0045 bool filter(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0046 const edm::EDGetTokenT<reco::VertexCollection> vertexSrc;
0047 const unsigned int minNDOF;
0048 const double maxAbsZ;
0049 const double maxd0;
0050
0051 };
0052
0053 GoodVertexFilter::GoodVertexFilter(const edm::ParameterSet& iConfig)
0054 : vertexSrc{consumes<reco::VertexCollection>(iConfig.getParameter<edm::InputTag>("vertexCollection"))},
0055 minNDOF{iConfig.getParameter<unsigned int>("minimumNDOF")},
0056 maxAbsZ{iConfig.getParameter<double>("maxAbsZ")},
0057 maxd0{iConfig.getParameter<double>("maxd0")} {}
0058
0059 GoodVertexFilter::~GoodVertexFilter() {}
0060
0061 bool GoodVertexFilter::filter(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0062 bool result = false;
0063 edm::Handle<reco::VertexCollection> pvHandle;
0064 iEvent.getByToken(vertexSrc, pvHandle);
0065 const reco::VertexCollection& vertices = *pvHandle.product();
0066 for (reco::VertexCollection::const_iterator it = vertices.begin(); it != vertices.end(); ++it) {
0067 if (it->ndof() > minNDOF && ((maxAbsZ <= 0) || fabs(it->z()) <= maxAbsZ) &&
0068 ((maxd0 <= 0) || fabs(it->position().rho()) <= maxd0))
0069 result = true;
0070 }
0071
0072 return result;
0073 }
0074
0075
0076 DEFINE_FWK_MODULE(GoodVertexFilter);