File indexing completed on 2024-04-06 12:33:18
0001
0002
0003
0004 #include "FWCore/Framework/interface/Event.h"
0005 #include "FWCore/Framework/interface/MakerMacros.h"
0006 #include "FWCore/Framework/interface/global/EDProducer.h"
0007 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0008 #include "FWCore/Utilities/interface/InputTag.h"
0009
0010 #include "CommonTools/Utils/interface/StringCutObjectSelector.h"
0011
0012 #include "DataFormats/TrackReco/interface/Track.h"
0013 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0014 #include "DataFormats/VertexReco/interface/Vertex.h"
0015 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0016
0017 #include <algorithm>
0018 #include <numeric>
0019 #include <memory>
0020 #include <vector>
0021
0022
0023
0024
0025 class TrackFromPVSelector : public edm::global::EDProducer<> {
0026 public:
0027 explicit TrackFromPVSelector(edm::ParameterSet const& iConfig);
0028
0029 void produce(edm::StreamID, edm::Event& iEvent, edm::EventSetup const& iSetup) const override;
0030
0031 private:
0032 double max_dxy_;
0033 double max_dz_;
0034 edm::EDGetTokenT<std::vector<reco::Vertex>> v_recoVertexToken_;
0035 edm::EDGetTokenT<std::vector<reco::Track>> v_recoTrackToken_;
0036 };
0037
0038
0039
0040
0041
0042
0043 TrackFromPVSelector::TrackFromPVSelector(edm::ParameterSet const& iConfig)
0044 : max_dxy_{iConfig.getParameter<double>("max_dxy")},
0045 max_dz_{iConfig.getParameter<double>("max_dz")},
0046 v_recoVertexToken_{consumes<std::vector<reco::Vertex>>(iConfig.getParameter<edm::InputTag>("srcVertex"))},
0047 v_recoTrackToken_{consumes<std::vector<reco::Track>>(iConfig.getParameter<edm::InputTag>("srcTrack"))} {
0048 produces<std::vector<reco::Track>>();
0049 }
0050
0051
0052
0053
0054
0055
0056 void TrackFromPVSelector::produce(edm::StreamID, edm::Event& iEvent, edm::EventSetup const&) const {
0057 edm::Handle<std::vector<reco::Vertex>> vertices;
0058 iEvent.getByToken(v_recoVertexToken_, vertices);
0059
0060 edm::Handle<std::vector<reco::Track>> tracks;
0061 iEvent.getByToken(v_recoTrackToken_, tracks);
0062
0063 auto goodTracks = std::make_unique<std::vector<reco::Track>>();
0064 if (!vertices->empty() && !tracks->empty()) {
0065 auto const& vtxPos = vertices->front().position();
0066 std::copy_if(
0067 std::cbegin(*tracks), std::cend(*tracks), std::back_inserter(*goodTracks), [this, &vtxPos](auto const& track) {
0068 return std::abs(track.dxy(vtxPos)) < max_dxy_ && std::abs(track.dz(vtxPos)) < max_dz_;
0069 });
0070 }
0071 iEvent.put(std::move(goodTracks));
0072 }
0073
0074 DEFINE_FWK_MODULE(TrackFromPVSelector);