File indexing completed on 2023-10-25 10:02:06
0001
0002
0003
0004
0005
0006
0007
0008 #include "FWCore/Framework/interface/stream/EDFilter.h"
0009 #include "FWCore/Framework/interface/EventSetup.h"
0010 #include "FWCore/Framework/interface/Event.h"
0011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0012 #include "FWCore/Framework/interface/ConsumesCollector.h"
0013
0014 #include <FWCore/ParameterSet/interface/ConfigurationDescriptions.h>
0015 #include <FWCore/ParameterSet/interface/ParameterSetDescription.h>
0016
0017 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0018 #include "DataFormats/VertexReco/interface/Vertex.h"
0019
0020 #include <algorithm>
0021
0022 class RecoTauPileUpVertexSelector : public edm::stream::EDFilter<> {
0023 public:
0024 explicit RecoTauPileUpVertexSelector(const edm::ParameterSet& pset);
0025 ~RecoTauPileUpVertexSelector() override {}
0026 bool filter(edm::Event& evt, const edm::EventSetup& es) override;
0027 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0028
0029 private:
0030 edm::InputTag src_;
0031 double minPt_;
0032 bool filter_;
0033 edm::EDGetTokenT<reco::VertexCollection> token;
0034 };
0035
0036 RecoTauPileUpVertexSelector::RecoTauPileUpVertexSelector(const edm::ParameterSet& pset)
0037 : minPt_(pset.getParameter<double>("minTrackSumPt")) {
0038 src_ = pset.getParameter<edm::InputTag>("src");
0039 token = consumes<reco::VertexCollection>(src_);
0040 filter_ = pset.getParameter<bool>("filter");
0041 produces<reco::VertexCollection>();
0042 }
0043
0044 bool RecoTauPileUpVertexSelector::filter(edm::Event& evt, const edm::EventSetup& es) {
0045 edm::Handle<reco::VertexCollection> vertices_;
0046 evt.getByToken(token, vertices_);
0047 auto output = std::make_unique<reco::VertexCollection>();
0048
0049 if (vertices_->size() > 1) {
0050
0051
0052
0053 std::remove_copy_if(vertices_->begin() + 1, vertices_->end(), std::back_inserter(*output), [this](auto const& vtx) {
0054 double trackPtSum = 0.;
0055 for (reco::Vertex::trackRef_iterator track = vtx.tracks_begin(); track != vtx.tracks_end(); ++track) {
0056 trackPtSum += (*track)->pt();
0057 }
0058 return trackPtSum > this->minPt_;
0059 });
0060 }
0061 size_t nPUVtx = output->size();
0062 evt.put(std::move(output));
0063
0064 if (!filter_)
0065 return true;
0066 else
0067 return nPUVtx;
0068 }
0069
0070 void RecoTauPileUpVertexSelector::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0071
0072 edm::ParameterSetDescription desc;
0073
0074 desc.add<double>("minTrackSumPt");
0075 desc.add<edm::InputTag>("src");
0076 desc.add<bool>("filter");
0077
0078 descriptions.add("recoTauPileUpVertexSelector", desc);
0079 }
0080
0081 #include "FWCore/Framework/interface/MakerMacros.h"
0082 DEFINE_FWK_MODULE(RecoTauPileUpVertexSelector);